UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

A 1D ocean mixing model of the Strait of Georgia : ecological responses to physical forcing 2006

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata


ubc_2006-0023.pdf [ 11.7MB ]
JSON: 1.0052593.json
JSON-LD: 1.0052593+ld.json
RDF/XML (Pretty): 1.0052593.xml
RDF/JSON: 1.0052593+rdf.json
Turtle: 1.0052593+rdf-turtle.txt
N-Triples: 1.0052593+rdf-ntriples.txt

Full Text

A ID Ocean Mixing Model of the Strait of Georgia Ecological Responses to Physical Forcing by A. Kathleen Collins . B.Sc, Mount Allison University, 2002 A THESIS SUBMITTED IN PARTIAL F U L F I L M E N T OF THE REQUIREMENTS FOR T H E D E G R E E OF M A S T E R OF SCIENCE in The Faculty of Graduate Studies (Oceanography) THE UNIVERSITY OF BRITISH COLUMBIA December, 2005 © A. Kathleen Collins, 2005 11 Abstract A coupled biophysical model of the Strait of Georgia, British Columbia, Canada has been devel- oped and successfully predicts the timing of the spring phytoplankton bloom. The physical model is a one dimensional vertical mixing model that is forced with hourly winds, cloud fraction, air temperature, and humidity, daily hydrographic data, and initial profiles of temperature, salinity, and fluorescence. The physical model uses a K Profile Parameterization of the boundary layer, and includes local parameterizations for albedo, cloud filtering, light attenuation, heat and fresh- water fluxes. The Fraser River is used to parameterize freshwater flux, horizontal advection, and upwelling. The biological model includes one phytoplankton class (microphytoplankton) and one nutrient source (nitrate). The coupled biophysical model was tested to determine what physical factors are controlling the arrival of the spring bloom. Wind was found to control the timing of the bloom arrival, with strong winds delaying the bloom and weak winds causing the bloom to arrive earlier. Solar irradiance had a small effect on the arrival time and freshwater input (primarily Fraser River discharge) was insignificant to the arrival time of the bloom. iii Table of Contents Abstract ii Table of Contents iii List of Tables vii List of Figures viii Acknowledgements xi 1 Introduction 1 1.1 Layout and Objective of Document 1 1.2 Physical Description of Strait of Georgia 2 1.2.1 Overview 2 1.2.2 Eraser River 3 1.2.3 Winds 7 1.2.4 Tides and Currents 8 1.3 Biological Description of Strait of Georgia 9 1.3.1 Bloom Basics 9 1.3.2 Phytoplankton Community 9 Table of Contents iv 1.3.3 Nutrients 11 1.3.4 Role of Zooplankton and Upper Trophic Levels 13 1.4 Literature Review 14 1.4.1 Previous studies related to S T R A T O G E M 14 1.4.2 Related Biophysical Studies 18 2 M e t h o d s 22 2.1 Field Data 22 2.2 The Physical Model , 24 2.2.1 Theory 24 2.2.2 External Forcing 30 2.2.3 Bottom Boundaries 50 2.3 The Biological Model 51 2.3.1 Theory 51 2.3.2 Bottom Boundaries 55 3 Results 56 3.1 Physical Model Results 56 3.1.1 Model Tuning and Sensitivities . 56 3.1.2 Vertical Variability 58 3.1.3 Temporal Variability 66 3.1.4- Testing Model Performance 75 3.2 Biological Model Results 79 3.2.1 Model Tuning and Sensitivities 79 Table of Contents v 3.2.2 Vertical Variability 81 3.2.3 Temporal Variability 82 3.2.4 Testing Model Performance 85 4 Discussion 93 4.1 Interannual Variability in the Strait of Georgia 93 4.1.1 Wind 93 4.1.2 Fraser 94 4.1.3 Clouds 97 4.2 What is controlling the spring bloom? 98 4.2.1 Winds 99 4.2.2 Clouds 100 4.2.3 Fraser River 103 4.3 Tides 105 4.4 Model Performance Summary I l l 4.4.1 Model Limitations I l l 4.4.2 Consistency of Results 114 4.4.3 Comparison with Other Results 115 4.5 Implications 117 4.5.1 Climate 117 4.6 Future Work 119 5 Conclusion 120 Bibliography 123 ) Table of Contents vi A Physical Model 131 B Biological Model 137 V l l List of Tables A . l Parameters used in the Physical Model 135 A. 2 New cloud model coefficients 136 B. l Parameters used in the Biological Model 138 B.2 Model runs used to predict the spring blooms of 2002, 2003, 2004, and 2005 139 B.3 Model tests identification scheme 139 V l l l List of Figures 1.1 Region of study shown on map of North America 3 1.2 Map of Strait of Georgia! 4 1.3 Contours of observed salinity and temperature 5 1.4 Fraser River flow 6 1.5 Wind speed and direction histogram from Jan 2001 to July 2005 7 1.6 Tides at Point Atkinson 8 1.7 Contours of observed chlorophyll and nitrate.' • • • 10 1.8 Fraser River water quality parameters. . ; 12 1.9 Seasonal variation of zooplankton • • 13 2.1 Location map of STRATOGEM field sampling program 23 2.2 Englishman River flow at Parksville and precipitation at Nanaimo, 2002. . . . . . . 32 2.3 Fraser and Englishman River flow; surface salinity fit 34 2.4 The pressure gradient scheme , 38 2.5 The light system in the Strait of Georgia 41 2.6 Experimental set-up for albedo calculation 42 2.7 Model calculated PAR and cloud fraction 44 2.8 Irradiance parameterization fit choices 47 List of Figures ix 2.9 Nutrients at S3 52 3.1 Profiles of physical model parameters. May 9- June 6, 2002 61 3.2 Profiles of physical model parameters. December 4 2002- January 10 2003 62 3.3 Profiles of physical model parameters. April 23- May 28, 2003 63 3.4 Profiles of physical model parameters. July 3- August 13, 2003 64 3.5 Profiles of physical model parameters. March 2- March 17, 2004 65 3.6 Eddy diffusivity, eddy viscosity and density. 66 3.7 Mixed layer depth, wind magnitude cubed, and Eraser River flow 68 3.8 Surface salinity and Eraser River flow 70 3.9 Contours of modelled temperature and salinity 71 3.10 Model computed surface temperature and air temperature 72 3.11 Modelled IPAR compared to observed PAR 73 3.12 Modelled temperature at 3 m and temperature measured by the ferry 77 3.13 Modelled 3 m salinity and salinity measured by the ferry 78 3.14 Contours of modelled chlorophyll and nitrate 83 3.15 Modelled and observed profiles of phytoplankton and nitrate 84 3.16 Modelled surface phytoplankton and bottle sampled chlorophyll 86 3.17 Modelled surface nitrate and bottle sampled nitrate 87 3.18 Effect of wind on chlorophyll throughout the water column . . . 89 3.19 Critical depth and phytoplankton 90 3.20 Modelled 3 m phytoplankton and ferry recorded chlorophyll 91 4.1 Wind magnitude cubed for the whole sampling period 94 List of Figures x 4.2 Wind magnitude cubed and modelled surface chlorophyll 95 4.3 Fraser River flow for 2001-2005 96 4.4 Fraser River flow and surface chlorophyll 96 4.5 Modelled daily maximum irradiance for 2002-2005 97 4.6 Irradiance in spring 2002 and 2005. . . 98 4.7 Effect of different wind magnitude on 2003 bloom timing 101 4.8 Effects of the wind magnitude on the timing of the spring bloom 102 4.9 Effect on bloom arrival of testing the model with different amounts of cloud cover for 2004.104 4.10 Comparison of true bloom arrival date to predicted 105 4.11 Effects of cloud cover on the timing of the spring bloom 106 4.12 Effect on bloom arrival of testing the model with different amounts of flow for 2003. 107 4.13 Effects of Fraser River flow on the timing of the spring bloom 108 4.14 Tidal excursion and ferry recorded salinity 110 A . l Upwelling mechanism in the Strait 134 A.2 Region of upwelling in the Strait 134 A.3 Cloud fits 136 xi Acknowledgements This thesis is a small part of a larger project, and the work I present here has been accomplished only with the help and support of many other people, to whom I express acknowledgements here. First I would like to thank my supervisor, Susan Allen, for her patience, insights, and all-around unending support. Susan has been an excellent mentor to me throughout my degree, and I have learned many things from her- most importantly that ocean modelling can be exciting. I would also like to thank Rich Pawlowicz for always asking the hardest questions, and Mark Halverson and Olivier Riche for their cooperation and information sharing. I'd also like to thank those three guys for all their hard work on the STRATOGEM Field Program, without which I would have no thesis. Thanks also to Akash Sastri, Rana El-Sabaawi, and Maite Maldonado for their support in the biological realm of things. I'd like to acknowledge the people who have lead me in this direction and supported me along the way. Thanks to Kate Moran for opening the door to Oceanography, and for keeping it open. For all their support I thank my parents- for calling every Sunday night and always critiquing my papers and convincing me, way back when, that Math was not the enemy. Thanks to Adam for joining me on the wetter coast and being both a wonderful distraction and a constant sounding board. Chapter 1. Introduction 2 The Introduction chapter of this thesis includes a Literature Review, a physical and a bio- logical description of the study area. Next, the Methods chapter describes the parameters of the model in use, such as how the physical and biological models were adapted to the Strait of Georgia, followed by a Results chapter which includes physical and biological results of the modelling. A Discussion chapter answers the question of what physical factors are important in the arrival of the bloom and possible implications of the results. The thesis finishes with a Conclusions chapter, and Appendices follow the main text. 1.2 Physical Description of Strait of Georgia 1.2.1 Overview The Strait of Georgia is a semi-enclosed estuary located off the coast of British Columbia, Canada (Figure 1.1), bounded by mainland BC on the eastern side and Vancouver Island on the western side (Figure 1.2). It is 2xl0 5 m long, ranging in width from 2xl0 4 - 4xl0 4 m, and depths up to 400 m (Foreman et al., 1995). The Strait is connected to the Pacific Ocean by the Juan de Fuca Strait via Haro Strait, a well-mixed, narrow and shallow passage (Figure 1.2). Circulation in the estuary can be considered to be two-layer flow, witli less dense water flowing out on the top layer and more dense water flowing in at depth, while temperature follows a typical seasonal pattern (Figure 1.3). An annual cycle of Deep Water Renewal brings nutrients into the Strait of Georgia in the fall. Chapter 1. Introduction 3 150°W 125°W 100°W 75°W 60°W Figure 1.1: Region of study shown on map of North America. 1.2.2 Fraser River Seventy percent of the freshwater entering the Strait comes from the Fraser River (Waldichuck, 1957). The Fraser River is the largest river in British Columbia, draining an area one third of the size of the province (Canadian Heritage Rivers System, 2005). From its headwaters on the western side of Mount Robson in the Rocky Mountains, the Fraser River travels 1.4xl06 m to the Strait of Georgia, carrying an annual silt load of approximately 20 million tonnes of silt, clay, and gravel, 16.5 million tonnes of which is deposited into the Strait of Georgia (Canadian Heritage Rivers System, 2005). Flow is low during the winter months (average winter flow is 8.5xl0 2 m 3 s _ 1 (BC Rivers, 2005)), begins to increase in the early spring and then peaks in June when the snow and glaciers in the mountains melt; this is termed the freshet (average peak flow is 8x l0 3 m 3 s _ 1 (BC Rivers, 2005)). Thus, the volume and timing of the peak flow depends on climate: temperature, sunlight, rainfall, and size of snow pack. Figure 1.4 shows the annual mean flow of the Fraser River measured at Hope, approximately 1x10s m upstream of where the Fraser River enters the Strait of Georgia as a riverine plume (Figure 1.2). Chapter 1. Introduction 4 Figure 1.2: Map of Strait of Georgia, Strait of Juan de Fuca, Haro Strait, Vancouver and Vic- toria, and locations of data acquisition: Point Atkinson (tide), Sand Heads (wind), Hope (Fraser River), Y V R (meteorological variables), Parksville (Englishman River), Nanaimo, and Halibut Buoy; location of modelled sampling station is mid-Strait of Georgia. Chapter 1. Introduction 5 5 _ 10 E B 15 o. 8 20 25 30 ^7 v May Sep 2003 May Sep 04 May Sep 05 J I I w If May Sep 2003 May Sep 04 May Sep 05 Figure 1.3: Contour of salinity (upper panel) and temperature (lower panel, °C) measured mid- Strait on S T R A T O G E M monthly cruises for April 2002 - June 2005. The fresh water signal of the Fraser River is apparent in the salinity plot, and the seasonal cycle of the temperature is evident in the temperature plot. Chapter 1. Introduction (i 16000F Day of Year Figure 1.4: Mean (thick line), maximum (upper thin line), and minimum (lower thin line), and one standard deviation on each side of the mean (grey shading) of annual Fraser River flow from 1912 to 2004. When the fresh water enters the Strait of Georgia it is known as the Fraser River plume, a large volume of river water mixing with seawater. The plume is visible to the naked eye because the high silt content of the Fraser River makes the plume appear as a large brown area beginning at the mouth of the Fraser River and extending up to ~1.5xl0 4 m North/South into the Strait of Georgia, and often across the Strait to the eastern shore of Vancouver Island. The edge of the plume is marked clearly by a front where the water changes from briny to saline. Mixing between seawater and the plume occurs when the river water entrains seawater from below. The large influx of potential energy associated with the plume stratifies and stabilizes the water column and causes the mixed layer depth to become shallower. Chapter 1. Introduction 7 Figure 1.5: Wind speed and direction from Jan 2001 to July 2005 measured at Sand Heads (left panel) showing the alignment of the strong winds with the axis of the Strait (right panel). The frequency of wind speed is depicted by the contours, so that it is clear that winds blow mainly from the NW and SE directions, and most of the time the wind is weak (<5ms _ 1). 1.2.3 Winds From late spring to early fall winds from the North and Northwest cause upwelling of high salinity water at the shelf near the western end of the Juan de Fuca Strait. From late fall to early spring winds from the Southeast favour coastal down-welling and low salinities. Strong winds in the Strait of Georgia blow along the axis of the Strait, but most of the wind is relatively weak (Figure 1.5) (Stronach, 1977). Wind plays an important role in the location and movement of the plume, especially at the Southern Arm of the Fraser River where the flow is not constricted by a jetty and can be blown across the Strait (Hoos and Packman, 1974). Wind can also break down the plume by causing mixing resulting in cold saline patches (Hoos and Packman, 1974; Halverson et al., 2005) . Chapter 1. Introduction 8 Figure 1.6: Upper panel: Tidal cycle for January 1 2003 to December 31 2003; lower panel: tidal cycle for January 1 2003 to February 1 2003 at Point Atkinson. 1.2.4 Tides and Currents Tidal currents are well studied in the Strait (LeBlond, 1983; Foreman et al., 1995). Tides are mixed semi-diurnal (Figure 1.6) and are forced by the Pacific Ocean tides. Tidal mixing is strongest in Haro Strait due to sills in the bathymetry and narrow constrictions of the islands. The tides advect the plume North on the flood tide and South on the ebb tide and also modulate the plume where it enters the Strait. Fraser River discharge is greatest at the large ebb, while the flood tide forms a salt wedge underneath the plume and hinders the discharge at the mouth of the Fraser River by pushing the river water back up the river (Stronach, 1977). When runoff conditions are low the wedge can intrude between 1.5xl04 and 3xl0 4 m inland (Swain, 1998). The movement of the plume is greatest when tides and wind act in concert. The surface water currents are mainly due to winds, tides, and the plume, with secondary movement accounted for by geostrophic forces such as the Coriolis force (Hoos and Packman, Chapter 1. Introduction 9 1974). The currents generated from winds travel at speeds of approximately 2% of the wind speed (Hoos and Packman, 1974). The Coriolis force causes the residual counter-clockwise circulation in the surface water (Stronach, 1977). 1.3 Biological Description of Strait of Georgia 1.3.1 B l o o m Basics A phytoplankton bloom typically occurs in the Strait of Georgia in the early spring, followed by an increase in zooplankton grazers and nutrient depletion in the summer (Figure 1.7) (Harrison et al., 1983). Phytoplankton growth is controlled by light in the winter and nutrient availability in the spring and summer. The size of the spring bloom is on average 30 mg m - 3 chlorophyll a, herein referred to as mg m - 3 chl. Populations below the surface layer are much smaller although of the same species as those in the surface layer (Cassis, 2005). Wind mixing in summer and fall can bring nutrients to the surface and cause smaller blooms to occur. 1.3.2 P h y t o p l a n k t o n C o m m u n i t y Species diversity in the Strait of Georgia follows a seasonal cycle with diversity at a minimum in the winter and at maximum during the summer (Hobson and McQuoid, 1997). The planktonic assemblage is composed mainly of flagellates in the winter and diatoms during the spring and summer (Cassis, 2005; Harrison et al., 1983), while during the spring bloom only a few species provide most of the cells. The spring bloom is often initially composed of the chain-forming diatoms Thalassiosira spy. and Skeletonema costatum and progresses with the explosion of S. costatum, while Thalassiosira spp. and the diatom Chaetoceros convolutus bloom more slowly and Chapter 1. Introduction 10 Figure 1.7: Contour of log chlorophyll (mgm 3 , upper panel) and nitrate (pM, lower panel) mea- sured to 30 m mid-Strait on STRATOGEM monthly cruises for April 2002 - June 2005. Surface nitrate fluctuations correspond to the bloom events, and high concentrations of nitrate are maintained at depth throughout the year. White patches on the nitrate plot indicate missing data due to incomplete nitrate samples for a given cruise. Chapter 1. Introduction 11 lag the 5. costatum bloom (Cassis, 2005; Hobson and McQuoid, 1997; Harrison et al., 1983). The progression from flagellates during the winter to diatoms during the bloom can be explained by the change in environmental conditions associated with the arrival of the bloom (Parsons et al., 1978) as well as by nutrient uptake kinetics (Parsons, 1979). Diatoms can out-compete flagellates in terms of maximum growth and dominate under spring bloom conditions of high nutrient availability and turbulence, while flagellates do better under low nutrient conditions and have lower maximum growth rates (Parsons et al., 1978; Parsons, 1979). Light limitation prevents the explosion of S. costatum during the winter. The influence of environmental conditions on the composition of the planktonic community was demonstrated by Price et al. (1985), who explained that differences in frontal and stratified phytoplankton communities in the Strait are due to differences in nutrient availability. Hobson and McQuoid (1997) observed that temperature, salinity, day length and irradiance, in decreasing order, were correlated to the temporal variation in cell species in the Strait, and that S. costatum was found throughout most of the year. 1.3.3 Nutrients Most of the nutrients in the Strait of Georgia enter from the Pacific Ocean through the Juan de Fuca Strait. Some of the nutrients are mixed into the upper layer by strong tidal mixing in Haro Strait, the rest continue with the estuarine flow into the Strait of Georgia. Wind and freshwater flux in the Strait of Georgia drive entrainment that brings the nutrients to the upper layer where the bloom occurs in the spring. Nitrate is abundant in the Strait throughout the water column in the winter, (Figure 1.7), and this nitrate is the primary source of nutrients to the phytoplankton during the spring bloom. Environmental conditions play a role in which type of nitrogen phytoplankton prefer (Price et al., Chapter 1. Introduction 12 12 F 1 - I Dissolved N03 (u M) Figure 1.8: Water properties for Fraser River measured at Hope during 1996, 1997, and 1998 (BC Rivers, 2005). Note: NTU=Nephelometric Turbidity Units. Left axis: dissolved NO3. Right axis: turbidity and flow. 1985). In the Strait of Georgia, nitrate is the primary source of nitrogen as levels of ammonium and urea are very low (Harrison et al., 1983). The possibility of eutrophication due to anthropogenic input of nitrates is low because a balance exists between nitrate entering the Strait through the Juan de Fuca Strait and nitrates leaving with the estuarine circulation (Mackas and Harrison, 1997). Nitrogen levels in the Fraser River are low (Figure 1.8); in fact, the open ocean source for nitrates is an order of magnitude larger than the total freshwater source, which includes the Fraser River (Mackas and Harrison, 1997). Eutrophication will not occur without a change in anthropogenic sources of nitrate or a change in the time of year when anthropogenic sources of nitrate enter the Strait (Mackas and Harrison, 1997). Chapter 1. Introduction 13 60 — Copepods dry wBight 0-100m (ug/L) >20um Chlorophyll, 0m (mg m" 3 ) 50 h 40 30 20 10 k " 1 Sep Jan May 2003 May Sep 04 May Sep 05 May Sep Figure 1.9: Seasonal variation of zooplankton in the Strait of Georgia. Copepods typically peak 1.3.4 R o l e of Zoop lank ton and U p p e r Trophic Levels The surface macrozooplankton community in the Strait of Georgia is mainly composed of copepods, and is dominated by the large calanoid copepod Neocalanus plumchrus (Figure 1.9) (Harrison et al., 1983). Copepod abundance typically peaks following the spring phytoplankton bloom. Euphausi- ids, chaetognaths, and some deep-living copepods live in mid-depths and deep water (Harrison et a l , 1983). The copepods E. bungii, Calanus marshellae, C. pacificus and Pseudocalanus minutus are also important species in the Strait at different times of the year (Harrison et al., 1983). The sur- vival and growth of N. plumchrus depends on phytoplankton biomass when the young copepods migrate to the surface in March or April of every year (Kobari and Ikeda, 2001). One species of primary importance to the commercial fishery is salmon, a fish that preys heavily on zooplankton such as copepods (Parsons et al., 1984). Other commercially and recreationally important fish in the Strait of Georgia include the lingcod, rockfish, and herring, in addition to up to seven salmon species (Fisheries and Oceans, Canada, 2005), all of which obtain their food following the phytoplankton peak. Chapter 1. Introduction 14 from some level of the food web of the Strait of Georgia. Exactly how much these fish depend on the lower trophic levels such as plankton, and therefore how much they depend on environmental changes is a focus of STRATOGEM. 1.4 Literature Review 1.4.1 Previous studies related to STRATOGEM The wheels of S T R A T O G E M were informally set in motion at the Strait of Georgia Workshop in November, 1972 at the Institute of Oceanography at UBC. Leading regional scientists met to discuss the physics, biology, and human usage of the Strait, and how best to study, model and understand the system in the future. One of the main goals presented at the meeting was to "measure certain time series which will clarify the level of natural variation through which we hope to perceive some trend," (Gargett, 1972), a goal echoed thirty years later by the S T R A T O G E M project. The conclusion of the meeting was that "we do not at present know enough about the statistics of any variable in order to suggest a sensible program for monitoring water quality in the Strait of Georgia," (Gargett, 1972). The report stated that the Strait was "full of problems which invite solutions," (Gargett, 1972). This section will review some of the attempted solutions and some of the steps that have been taken toward reaching a solution since the workshop. Physical Oceanography as a science grew enormously in the late 1970s, allowing for the first comprehensive report of the modern age of the Strait (LeBlond et al., 1983). New technology allowed for a higher frequency of observations, new data analysis techniques, and the result was a very in-depth look at the physical processes occurring in the Strait. LeBlond concluded his report by stating that a working upper layer model was the next step towards understanding the Chapter 1. Introduction 15 Strait, but that observational data at high and low frequency scales were needed to achieve this goal. LeBlond also focused attention on the importance of deep water renewal to the ecology of the Strait and the need for it to be well-modelled. Although his findings were highly significant LeBlond acknowledged that the new information in his report was only a fraction of the complete picture of the "inner workings" of the Strait, and that much more data, comprehensive and working models, and theoretical understanding were required in order to move forward. LeBlond's study set the basic framework for this thesis. One of the first studies to investigate biology and physics in the Strait of Georgia was the model by Stronach (1977), who attempted to relate biological properties of the Strait of Georgia to mixing and plume dynamics of the Eraser River using a model of Fraser River estuary. He reported velocity differences between the plume and the underlying water of up to 3 m s - 1 at the mouth, causing significant vertical mixing and introduction of nutrients into the upper layer (Stronach, 1977). His model showed that tidal cycles play a crucial role in understanding plume dynamics, but it lacked wind forcing and therefore left room for improvement in the modelling of biophysics of the Strait of Georgia. In his model relating biology and physics in the Strait, De Lange Boom (1976) showed that an increase in solar radiation, a decrease in the upper layer thickness due to greater freshwater input, and an increase in maximum production rate all increase production in spring but the simultaneous increase in turbidity and advection decreases available light in the water column and acts to decrease production. He concluded that the effects of physical processes on biological processes are greater than the converse- the advection of chlorophyll is more important than light absorption by plankton to the overall spatial and temporal distribution of production. De Lange Boom's model had a poor advection scheme, however, and neither did it represent short-term Chapter 1. Introduction 16 effects on the ecosystem such as diurnal variations in insolation, tidally induced variations in the velocity field, or migration of zooplankton. Parsons and Kessler (1987) used an ecosystem model to analyze the effects of the timing of a plankton bloom on young salmonoid growth and survival in an estuary similar to the Strait of Georgia. They assessed changes in mixed layer depth, temperature, silt load and thus light penetration depth, total surface radiation, changes in ctenophore biomass (competition for zoo- plankton prey), changes in initial biomass of salmon and salmon arrival lag times. They found that the physical environmental conditions and the initial standing stock.of ctenophores can both have a direct effect on the biomass of salmon through changes in the standing stock of zooplankton. In addition, they found that it is the quality of zooplankton growth conditions that determines whether or not initial biomass of salmon and timing of salmon arrival to a spring bloom already in progress have an effect on the ultimate growth of salmon (Parsons and Kessler, 1987). A large amount of work related directly to linking physics and the spring bloom was undertaken by Yin et al. in the 1990s (1997b; 1996). Through a series of cruises a data set containing biological and physical observations was obtained for the Strait of Georgia and the role of wind, Fraser River discharge, and zooplankton grazing on the spring phytoplankton dynamics were investigated. Yin et al. (1997b) stated that a stratified water column is necessary for a bloom to occur, and that this stratification is achieved by freshwater input and destroyed by winds greater than 4ms - 1 . They further asserted that the bloom will, in fact, not fully develop until the freshwater discharge is strong enough to overcome wind and tidal mixing (Yin et al., 1997b) and that the bloom arrival is dependent on the freshet arrival (Yin et al., 1996). Y in (1997b) concluded that "interannual variations in winds and the timing of the annual freshet determine the timing and duration of the spring bloom". Furthermore, Yin (1997b) asserted that after the initial bloom has begun, wind Chapter 1. Introduction 17 can interrupt the bloom, and when the wind event is over the surface N O 3 will have increased and the bloom will continue; this same mechanism of nitrate upwelling and subsequent chl a increase is observed during the summer when surface concentrations are low (Yin et al., 1997a). A two-part study on seasonal and interannual variability relating estuarine circulation to food web dynamics in the Strait was undertaken.by L i et al. (1999; 2000). A two-layer box model of the physical system was developed (Li et al., 1999) connecting the Strait of Georgia to the Juan de Fuca Strait through Haro Strait. The first part of the study investigated the physics of the system. It used a model forced primarily by a hyperbolic function representing Fraser River runoff and also included vertical mixing between layers and the spring-neap cycle of tidal mixing between basins. In the second part of the study the box model was coupled to an ecosystem model in which concentrations of nutrients and populations of phytoplankton and zooplankton were varied over seasonal and interannual cycles (Li et al., 2000). Phytoplankton growth in the model was controlled by one of nutrient or light limitation, both of which are governed by Michaelis-Menton equations, and phytoplankton biomass is controlled by zooplankton grazing. Mortality rates of plankton were set to be proportional to population size and phytoplankton and nutrients were allowed to move vertically and horizontally through boxes while zooplankton were assumed to be able to keep themselves in the sunlit upper layer and therefore were confined to the upper boxes (Li et al., 2000). The studies were an attempt to investigate possible links between physical regime variations to changes in the planktonic ecosystem, and subsequently to the upper trophic levels and the fisheries. The modelled estuarine circulation responded quickly and in accordance with observations to decadal changes in Fraser River runoff volume associated with climatic regime shifts related to mountain snow accumulation. Plankton populations, however, were found to be insensitive to Chapter 1. Introduction 18 these interannual changes in the estuarine circulation; L i et al. (2000) suggest that plankton might actually respond to climate variability through changes in their biological rate parameters . This thesis will continue the investigation into what factors cause variability in plankton populations. 1.4.2 Re la t ed B iophys i ca l Studies The number of studies attempting to link physical forcings to ecological responses has been in- creasing over the past decade. Some of the links that have been discovered put S T R A T O G E M at the forefront of the research community; a few of these studies are discussed below. Beamish et al. (1997) showed that abundance of Fraser River sockeye salmon is related to changes in climate trends of the North Pacific ocean-atmosphere system . Edwards and Richardson (2004) and Richardson and Shoeman (2004) strengthened the link between climate change and the lower and upper trophic levels by showing that "not only is the marine pelagic community responding to climate changes, but also that the level of response differs throughout the commu- nity and the seasonal cycle, leading to a mismatch between trophic levels and functional groups" (Edwards and Richardson, 2004). This change in the synchrony of timing between primary, sec- ondary, and tertiary production is already having an effect on the marine trophodynamics and on commercially important fish species in the North Atlantic; one example is the decline in the North Sea cod stocks (Edwards and Richardson, 2004; Richardson and Shoeman, 2004). There is strong evidence for a link between the increasing trend in Northern Hemisphere temperature and North Atlantic Oscillation and the bio-geographical shifts of North Atlantic copepod species (Beaugrand et al., 2002). It is links such as these that S T R A T O G E M ultimately seeks to uncover, the first step being to elucidate the effect of climate change on the primary production. The international and inter-oceanic study GLOBEC, Global Ocean Ecosystem Dynamics, was Chapter 1. Introduction 19 formed in 1991 as a multi-disciplinary research program designed by oceanographers, fisheries scientists, and marine ecologists to examine the potential impact of global climate change on ocean ecosystems (Batchelder, 2003). The study, which is partially complete but still ongoing, specifically aimed to investigate the distribution, abundance, and production of marine animals, from plankton to commercial fish species, using diagnostic and prognostic ecosystem models capable of capturing the ecosystem response to major climatic fluctuations (Batchelder, 2003). Jin et al. (2005) developed a coupled biophysical model of the Bering Sea middle shelf and used it to correctly predict the spring phytoplankton bloom and investigate effects of physical forcings on the bloom dynamics. The model is forced with cloud cover, precipitation rate, sea level pressure, relative humidity, wind speed, air and sea temperature, and salinity; fluorescence measurements were used for model comparison. Tidal mixing, wind stirring, and stratification were tested for their effect on the timing, duration and magnitude of the bloom. Absence of tidal mixing caused a decreased, shorter bloom, whereas absence of wind and vertical stratification caused an increase in mean phytoplankton concentration and longer blooms (Jin et al., 2005). Absence of tidal mixing did not affect the bloom timing whereas absence of wind and stratification caused the bloom to be delayed. Chen (2004) coupled a ID Nutrients, Phytoplankton, Zooplankton and Detritus (NPZD) model with a Finite Volume Coastal Ocean Model (FVCOM) of George's Bank to research the physical forcings of the ecosystem in both a shallow and a deep (>60 m) region. He was able to predict the spring bloom in the shallow, well-mixed region and found that the timing of the bloom depends on light intensity and attenuation. However, in the deep region he was not able to reproduce the spring bloom as successfully, but found that its timing is linked to the seasonal development of stratification. Chapter 1. Introduction 20 The spring blooms of 1993, 1995 and 1996 of Prince William Sound, Alaska were modelled successfully by Eslinger et al. (2001) using a ID mixed-layer simulation model coupled with an eight component biological model. The physical model is forced with local meteorological forcing from a moored buoy and uses observations of chlorophyll and nitrate from a field sampling program for biological comparison. The biological model uses a temperature dependent phytoplankton growth rate, and includes both diatoms and flagellates. Analysis of meteorological data and model simulated blooms led the authors to conclude that variations in bloom timing, rate of increase, and duration were due to differences in wind speed and air temperature in the early stages of stratification (Eslinger et al., 2001). Analysis of spring bloom dynamics in an estuary off the United Kingdom coast reveals that the main factor controlling the timing of the bloom is irradiance (Iriarte and Purdie, 2004). Linear regressions were conducted between chlorophyll a concentration and each of mean water column P A R (average of 7 days prior to sampling date), river flow, tidal range, and wind speed (daily mean averaged for 3 days prior to the sampling date) for 1988, 1992, 2001, 2002 and 2003. The topic of biological response to physical environmental forcings, therefore, has been well- studied with varying degrees of complexity and success. Although the objectives of the aforemen- tioned studies are similar to one another and to the objectives of S T R A T O G E M , the methodologies differ greatly. In the chapters that follow, I will develop a ID vertical biophysical model that pre- dicts the spring bloom for four very different - both physically and biologically- years. I will then use the model to detect any possible links between physical forcings and the timing of the spring bloom. Some of my methods will match those of the precursors to this study, but in many ways, my methodology will diverge from that of my contemporaries. The model will contain many local parameterizations, the likes of which have not been attempted by previous studies; these will be Chapter 1. Introduction 21 discussed in detail. I will use continuous data forcing and initialize the model with real data. I will tune the model to these four different years using only one unknown parameter, and I will use the model as a tool to investigate the interaction of physical forcings with the ecosystem. In these ways, this study will be significantly different from the studies reviewed above, and will provide new insight into the physical controls of the biological community of the Strait of Georgia that will hopefully answer some of the questions posed by Gargett (1972) more than thirty years ago. 22 Chapter 2 Methods 2.1 Field Data The field data used in this thesis were collected as part of the S T R A T O G E M project over four years at nine stations in the Strait of Georgia (Figure 2.1). The stations were selected so that data was gathered both inside and outside the the Fraser River plume and in the middle of the Strait as well as the Southern part. Sampling took place approximately once per month, and once per week during the spring bloom. Casts were taken with a Conductivity Temperature Depth (CTD) sensor; variables recorded at each station included fluorescence, Photosynthetically Active Radiation (Ip/t/j), temperature and salinity. Nitrate, phosphate, silicon, chlorophyll (for size fractionation analysis), salinity, oxygen, and taxonomy samples were taken at 0, 5, 10 and 30 m depth at each station, with additional deeper depths at some stations. Other field data include meteorological data obtained from Environment Canada's weather station at the Vancouver International Airport and at Sand Heads Meteorological Station (YVR, 2005). These include humidity, air temperature, cloud fraction data, and wind speed and direction. Fraser River and Englishman River flow data were obtained from the Water Survey of Canada (BC Rivers, 2005). Downward short-wave radiation data measured at the U B C Department of Agroecology was used for the parameterization of the cloud model that will be described in this Chapter 2. Methods 23 Figure 2.1: Location map of STRATOGEM field sampling program showing all nine sampling sta- tions and the ferry tracks (dashed lines). The central route (Tsawwassen to Departure Bay) crosses S3 and is the source of data used for comparing to model results. section (Ketler, 2005). Irradiance data from Halibut Bank buoy (Figure 1.2) was used to verify the cloud model (Gower, 2005). Part of the S T R A T O G E M project involves an innovative sampling program aboard several BC Ferries vessels. Two routes are traversed eight times daily, resulting in a comprehensive data set. The central route between Tsawwassen and Departure Bay (Figure 2.1) passes over Station S3 and measures temperature, salinity, fluorescence, latitude, longitude, and GPS time. Thus a large quantity of spatial and temporal data is available on a continuous timescale to compare with model output. Chapter 2. Methods 24 2.2 The Physical M o d e l 2.2.1 T h e o r y A one-dimensional vertical mixing model is used to model physical properties of the water column at S T R A T O G E M Station S3 (Figure 2.1, Latitude = 49 07.5 N , Longitude = 123 33.5 W) between 0 and 40 m depth. The model used is an adapted version of the K Profile Parameterization (KPP) nonlocal boundary layer model (Large et al., 1994). The model represents mixing due to turbulent vertical velocities of unresolved eddies and calculates numerical values for various physical properties in time and space. Briefly, the model theory is as follows: first the external forcing is prescribed, then diffusivity and nonlocal transport are computed and lastly the mixed layer depth is computed. In this section I will discuss each of these three steps beginning with the computation of diffusivity, then mixed layer depth. These are followed by a more detailed discussion of external forcings. A p p e n d i x A gives values for model parameters and specifics such as grid information. It is useful at this point to define and differentiate mixed layer depth (MLD), mixing layer depth, and the boundary layer. Mixing layer depth is defined by Thomson and Fine (2003) as the depth determined by a "balance between destabilizing effects of mechanical mixing and stabilizing effects of surface buoyancy flux". It is characterized by "uniform to near-uniform density, active vertical mixing, and high dissipation", and caused by turbulence generated by external forcing (Thomson and Fine, 2003). The mixed layer depth is of identical density to the mixing layer depth, but represents a time-integrated response to the previous mixing events and as such is usually deeper than the mixing layer. After the mixing stops, the turbulent eddies retreat to shallower depths; however, Thomson and Fine suggest that the erosion of stratification takes time, Chapter 2. Methods 25 and by the time the stratification has been eroded and a new mixed layer depth established, the eddies may have retreated so that the mixing layer depth is now above the mixed layer depth. The boundary layer referred to by Large et al in the K P P model is equivalent to mixed layer depth, not mixing layer depth. The model is initialized with profiles of CTD temperature (T), salinity (S), and chlorophyll fluorescence. Turbulent and non-turbulent surface fluxes are the forcing mechanisms that control mixing at the surface. Turbulent forcing includes wind stress components TU, TV, heat flux Qt (includes net long wave radiation, latent and sensible heat fluxes) and freshwater flux Ft. Non- turbulent forcing is derived from solar irradiance distributed through the water column Qn. These are described specifically in the External Forcing section. The external forcing is applied to the model continually as hourly meteorological data and daily river flow data. The external forcing on the surface of the water (the top of the model grid) introduces kinetic energy into the water column and causes mixing and the generation of turbulent eddies. The time evolution of a property X due to the eddies is expressed as the vertical divergence of the kinematic turbulent fluxes where X represents the property being computed (velocity, salinity, temperature, . . . ) , d refers to a partial derivative, t represents time, z is the upward vertical coordinate, x represents the fluctuations of a property X, w is the turbulent vertical velocity, and an overbar represents a time average. The depth dependent flux of turbulent eddies is determined by where "fx is a nonlocal transport term described below, and d is depth. Kx is the profile of boundary dtX = —dzwx (2.1) wx(d) = -Kx(dzX - 7X) (2.2) Chapter 2. Methods 26 layer diffusivity, prescribed by KxJa) = hwx{o)G{cj), (2.3) where h is mixed layer depth, wx is a generalized depth dependent turbulent velocity scale, a = d/h is a dimensionless vertical coordinate that varies from 0 to 1 in the mixed layer, and G(cr) is a non-dimensional vertical shape function. The nonlocal transport term jx accounts for the unstable (convective) situations when the length scale of the eddies may be equal to or greater than the boundary layer depth, and the local diffusion approach is no longer appropriate; diffusivity.is then better represented by a nonlocal diffusion term. Equation (2.2) dictates that the vertical turbulent flux is equal to the diffusion of the local (X) and nonlocal (7^) property gradients. The shape function G(a) is a cubic polynomial whose coefficients depend mainly on the mixed layer depth; scalar properties will have the same G(o) as each other but different from the non- scalar properties. The turbulent velocity scale wx depends on whether conditions in the bound- ary layer are convectively stable or unstable; a fixed profile is chosen based on these conditions (Large et al., 1994, p. 370). The mixed layer depth from the previous time step (h) is then used and the boundary layer diffusivity is calculated (2.3). From this the vertical turbulent fluxes of momentum and scalars are calculated (2.2). The internal mixing is calculated next. The ocean interior is defined as the section of the water column below h. Internal mixing is included in the model as shear instability mixing (v%(<£)) nor- malized by its maximum (vo), unresolved internal wave breaking shear (vx), and double diffusion due to salt fingering {vx{d)). The shear instability mixing calculation for vsx{d) is based on a comparison with the local Chapter 2. Methods 27 gradient Richardson number Rig = N2 (2.4) {dzUf + {dzVf where iV 2 is the local buoyancy frequency, U and V are mean velocities evaluated over all depths. A value for vx(d) is determined from a comparison of Rig with a critical value (equations 28a-c) (Large et al., 1994), and the result is a depth profile of shear instability mixing. Internal mixing due to double diffusion is calculated based on the double diffusion density ratio where a and (3 are the thermodynamic expansion coefficients for temperature and salinity, re- spectively (Large et al., 1994, equations 31a-31c), leading to a depth-dependent profile of double- diffusive mixing due to salt fingering, vx(d). Scalar and momentum quantities are not calculated separately for shear and double diffusion mixing. Following Large et al. (1994), both scalar and momentum quantities of mixing due to internal waves vx are assumed to be small and are treated as background diffusivities for the ocean interior. The profile of overall interior diffusivity is calculated by The values for the scalar and momentum quantities of v™ were chosen to optimize model output as described in the M o d e l Tuning and Sensitivities section of the Results chapter and are given in Table A . l . Two profiles now exist for the vertical turbulent fluxes- one in the mixed layer (Kx) and one in the interior (vx)\ these are matched together to create a continuous shape profile for the entire water column that is used in the semi-implicit integration to solve for new values of the properties, Rp = adzT/(3dzS (2.5) vx(d)=vx(d) + vxJ+vdx(d) (2.6) Chapter 2. Methods 28 X (Large et al., 1994, Appendix D). The turbulent vertical fluxes of momentum and scalars in the interior are then computed by wx(d) = -vx(d)dzX; ' (2.7) for all depths above and below the mixed layer depth (2.2 and 2.7, respectively) and the properties can be integrated in time and the depth profiles updated: dtU = dz[-wu} + fV-Py (2.8) dtV = dz[-wv]-fU-Px (2.9) dtT = dz[-^e}-dzQn (2.10) dtS = dz[-ws] - dzFn (2.11) where S is salinity, T is temperature, U and V are horizontal velocities, and / is the Coriolis pa- rameter, a constant. The parameterization for non-turbulent heat profile Qn is the solar irradiance profile (total light) discussed in the Light below the water section. The parameterization of the pressure terms px, py, which are not included in the open ocean scenario described by Large et al. (1994), is discussed in the Baroclinic Pressure Gradients section. The surface fluxes WUQ, WVQ, W6Q, and WSQ are computed following wuo = — - (2-12) Po l~v WVQ = (2.13) Chapter 2. Methods 29 % - (2.14) Po{0) where TU and T v . are the wind stresses are computed following (A.2) and are discussed in the W i n d section; So, po, and Cpo are the salinity, density, and specific heat at constant pressure of the surface seawater, respectively. The density of surface water of salinity 0 is denoted by po(0); Qt and Ft are the turbulent heat and freshwater fluxes, the parameterization of which are discussed in sections Light and Freshwater F lux sections, respectively. Calculations up to this point are completed with h at time t. The calculation of h at time t + 1 is based on the surface forcing, and on the oceanic buoyancy B(d) (2.16) and boundary layer horizontal velocity V(d) profiles B(d) = g(aT - j3S) (2.16) where a and j3 are evaluated at local values of temperature, T, and salinity, S. A near-surface reference buoyancy Br B^A lwJ°- {w^}  (2- 17) is calculated, where the subscripts 0 and hj, specify evaluation at the surface and at the bottom of the boundary layer, and I is total solar irradiance. A near-surface reference horizontal velocity vector Vr is also calculated. Br and Vr are estimates of the average buoyancy and velocity over the surface layer. The bulk Richardson number (Br-B(d))d Chapter 2. Methods 30 is calculated at every time step based on buoyancy Br (2.17) and velocity Vr of boundary layer eddies, and local buoyancy and velocity, B(d) (2.16) and V(d). Turbulent velocity shear is included as Vt/d. The mixed layer depth, h, is defined as the depth where Rib(d) (2.18) is less than the critical bulk Richardson number Ric, such that eddies penetrate to the depth where they become stable relative to local buoyancy and velocity. Ric is chosen to be 0.3 (Large et al., 1994). 2.2.2 External Forcing K P P (Large et al., 1994) was designed to model open ocean systems that differ from coastal systems such as the Strait of Georgia in several aspects. Many adaptations were made to the original model so that an estuarine system was being modelled and not an open ocean. These include adaptations to the horizontal processes governing the estuarine flow and depth-dependent processes governing light attenuation. Many of the models implicit to the K P P model were re-fit using local data in order to maximize the success of the model. The resulting model was optimized for the Strait of Georgia in order to be sensitive to the specific physical dynamics of the region. Firstly, the parameterization of freshwater was addressed. Large et al. (1994) include freshwa- ter input only as evaporation and precipitation. The importance of the Fraser River to the Strait of Georgia estuary requires a new riverine-based parameterization for freshwater flux. An estuarine entrainment process that depends on the influence of the Fraser River must also be modelled here, as it is absent from Large et al. (1994). Secondly, pressure gradients must be included as part of the estuarine flow. In the open ocean, flow in a given direction in a given layer can continue in a single direction without particular consequences. In a Strait it causes an accumulation of that layer downstream and a corresponding baroclinic pressure gradient against the flow; therefore Chapter 2. Methods 31 this process must be parameterized. Finally, the open ocean light absorption scheme of K P P was re-parameterized to model the coastal light system of the Strait of Georgia; this includes light in the atmosphere as well as light in the water. Many components of the model depend on a good model estimation of light- non-turbulent heat profiles are determined from the total light, and phytoplankton growth depends on the biologically available part of the light spectrum. These parameterizations are discussed in detail in the following sections. F r e s h w a t e r F l u x Physical dynamics in the Strait of Georgia are controlled by the freshwater input of the Fraser River, which drives the estuarine circulation. Daily Fraser River flow data measured at Hope, lx lO 5 m upstream of where the river meets the Strait is used as model input. In K P P freshwater flux is determined solely from evaporation and precipitation. However, in the Strait of Georgia the freshwater input from the Fraser River dominates the rain water input by many orders of magnitude. The inclusion of a riverine-based freshwater flux was one of the most important changes made to K P P to correctly model the Strait of Georgia. To accomplish this freshwater flux is calculated based on the change in surface salinity between timesteps, using a predicted salinity based on river flow data for the new salinity, and salinity from the previous model time step for the old salinity; the details of this process are as follows. Although the Fraser River flow dominates the freshwater flux into the Strait, effects of local precipitation must be included. Large rainfall events in the Georgia Basin are not always as severe upstream or at Hope (Figure 2.3) and the ferry data shows that significant rainfall events have a profound effect on the surface salinity in the Strait (Halverson et al.,- 2005). Thus, a local river, the Englishman River, was also used as part of the freshwater flux parameterization to Chapter 2. Methods 32 Figure 2.2: Englishman River flow at Parksville and precipitation at Nanaimo for 2002. The river flow is affected greatly by rainfall; flow and precipitation are highly correlated, with a coefficient of ~0.6. represent all other sources of freshwater input into the Strait. The Englishman River has its headwaters in the mountains of Vancouver Island, and due to the smaller snow-pack of these mountains it peaks in winter coincident with heavy precipitation and is lowest during the dry summer months. The river drains an area of 3.24xl08 m 2 over 3.8xl0 4 m directly into the Strait of Georgia at Parksville (124 16.58 W, 49 19.00 N) (Figure 1.2). Water flow of the Englishman River is greatly affected by rainfall and is subject to rapid fluctuation over short time periods (Mount Arrowsmith Biosphere Reserve, 2001). Although it is dammed the water is only released during the summer and early fall to supplement local water supply (City of Parksville, 2005), times when Fraser River flow is high and Englishman River flow is negligible in comparison. Englishman River flow correlated well with precipitation records at Nanaimo (Figure 2.2), so the Englishman River was the appropriate choice to include in the parameterization for freshwater flux. The first part of the freshwater flux parameterization was to fit the surface salinity at S3 to the Chapter 2. Methods 33 two sets of river flow data using a 3-day back average for the Fraser River (Qp) and a 19-day back average for the Englishman River (QE)- The back-averages are based on how long the freshwater signal stays in the water column and were chosen based on salinity records compared to storm records. For example, a high rain event on October 23, 2003 was still present at the surface on the November 6 cruise as a low peak of salinity (Figure 2.3), but Englishman River flow remained high only for a few days (Figure 2.2, 2.3). Model surface salinity calculations based on daily flow would therefore be too saline compared to observed salinity, which remained fresh for several weeks following the storm. A 19-day average of the Englishman River flow suitably modelled the surface salinity during periods of high flow and also accounted for the residence time of freshwater in the Strait. The length of time that the freshwater signal remains in the Strait is due to the staggered effect of the entrance of the storm water into the Strait. For example, an undammed river would flush immediately into the Strait but a dammed river or river that drains into a lake before the Strait would take several days to finally pour into the Strait. The dammed Coquitlam River is one such example; high runoff continues until Oct.31, 2003 following the aforementioned storm event. The Fraser River is much less variable than the local river system and a three day average was sufficient to model surface salinity successfully. The fit is Sriv = So — apQp — CLEQE (2-19) where the salinity at depth S£>=29.1166 and the fit coefficients are ai?=0.0019 s m - 3 , and a£=0.0392 The fit has R 2 = 0.67 (Figure 2.3b) and suitably models the surface salinity year-round and during storm conditions (Figure 2.3a). River data is only provided on a daily basis, giving one value per day of Sriv and leading to unstable model results. This requires that either Qp and Q E be interpolated in time to obtain a Chapter 2. Methods 34 Figure 2.3: (a) Fraser River and Englishman River flow for Jan. 2002 - June, 2005. (b) Averaged Fraser River and Englishman River flow data fit to surface salinity at S3. Arrows mark notable storm events. Chapter 2. Methods 35 river flow value every timestep, or that Sriv be smoothed over time. The latter option is chosen using a time constant of half a day, and the result is a predicted surface salinity based on river flow for each timestep: S(t + 1). The objective here is to obtain a value for how much freshwater was added to S3 between timesteps, and this can be calculated at this point as both the salinity at the current timestep, S(t + 1) and the modelled salinity from the previous timestep, S(t), are known. If we assume the surface layer of the model to be a box with height dz, we know that this box has total salt S(t)dz at time, t. This is equal to the total salt at the new salinity S(t + 1) of the original box dz plus the freshwater added, Ftdt S{t)dz = S{t + l){dz + Ftdt). (2.20) Equation 2.20 is solved for Ft, the freshwater flux, which is used to calculate the salinity flux (2.15). S(t)-S(t + l)dz Ft~ s(t + i) W  (2- 21) V e r t i c a l A d v e c t i o n o r U p w e l l i n g The estuarine circulation driven entrainment in the Strait dominates the vertical advection. As this vertical advection brings nutrients toward the surface we will refer to it as upwelling. Upwelling in the model is dependent on the forcing of the Fraser River, and upwelling velocity (w) is determined by calculating the influence of the Fraser River between S2-2 and S3 over the area between the stations (Figure A.2). It is assumed that the Strait is divided into two layers- a surface layer driven by Fraser River flow Qp and a deep layer with salinity Sp and flow QQ\ at S2-2 and QDI at S3, and similarly for surface parameters (Figure A . l ) . Entrainment (E) of water parcels from the lower layer into the upper layer is parameterized by Fraser River induced mixing. Wind is Chapter 2. Methods 36 not considered to contribute to entrainment of particles from the deep layer to the surface layer, because changes in mixed layer depth caused by wind mixing are calculated elsewhere in the model and are based on diffusivity rather than the process outlined here. The process of vertical advection parameterized here is driven entirely by estuarine circulation driven by the Fraser River. The Fraser River induced mixing is assumed to be the only entrainment between S3 and S2-2 and it occurs over some fraction 1/Cvv of the area between two concentric circles, Tr(r2i—r^2)/Cw, where r2 and r\ are the distances of each station to the mouth of the Fraser River (Figure A.2). Appendix A details the algebraic procedure that solves for E, and the end result is RSD(S2 - Si) (sD-s2)(sD-s1)  { • ' Similarly to the theory of (2.19), we know that the surface salinities are equal to the salinity at depth minus a factor proportional to the Fraser River flow, such that S\ = Sp — and S2 = So — Substituting, we get E = S o (72 — 7 1 ) , a n d upwelling velocity is defined entrainment over an area, so upwelling is calculated by w = CwSDi-n-rn), ( 2 . 2 3 ) ir(r%-r{) The value w is assumed to be the maximum upwelling velocity at 40 m depth. This value decreases linearly from the bottom to the surface giving a value for upwelling velocity for each grid level; the contribution of upwelled salt to the salinity profile is calculated following Jeffery (2002) (A.20) as detailed Appendix A. The upwelling of temperature, nitrate, and phytoplankton are calculated in the same way such that the Fraser River essentially controls all upwelling in the model. In Jeffery's model the upwelling system results in an imbalance of fresh and salt water at the surface as the only source of freshwater to the sub-Arctic Pacific is precipitation. She therefore prescribes a horizontal advection of freshwater into the surface to account for this problem; this Chapter 2. Methods 37 advection is not needed in the present model because the balance of fresh and salt water is achieved by the freshwater input of the Fraser River. Therefore there is no explicit horizontal advection in this model, and movement of water parcels is controlled by upwelling and diffusion. This procedure assumes that the upper layer flows, Q\ and Q2, exist in the upper 40 m and have salinities S\ and Perhaps a more appropriate value to represent upper layer salinity would be a value averaged over 40 m. The reader will find that using a rough estimate of this average, (Si + So)/2 results in the same estimate of E (2.22) given above, and so the assumption of surface salinity as upper layer salinity is not altogether unsuitable. The linear decrease of the upwelling velocity toward the surface may also strike the reader as overly simplistic. However, this solution produced results that were more favourable than results produced employing a more realistic scheme, such as a parabolic decrease of upwelling velocity. Future modelling efforts will continue to address this issue and attempt to find a more appropriate method to decrease upwelling velocity over depth. Baroclinic Pressure Gradients The open ocean is assumed horizontally large and so Large et al. (1994) do not include pressure terms in their model equations; but the enclosed nature of the Strait of Georgia requires their inclusion. To do this, it is assumed that (1) the basin is shaped like an ellipse with nearly vertical sides; (2) S3 is located in the exact center of the basin; (3) the basin is uniformly stratified (Figure 2.4). A flow in a given isopycnic layer will cause a thickening/thinning of that layer downstream / upstream (2.24) Chapter 2. Methods 38 S3 Om z ( m ) ci P. 40m 40m Figure 2.4: Cartoon depicting the scheme used to solve for the x and y baroclinic pressure gradi- where Ui is the baroclinic velocity Ui- < u'i >, Lx is the width of the ellipse, A z is the thickness between the two grid points and zxci is the thickness at the positive z-boundary of the ellipse. The thickness at the positive y-boundary is found similarly. These tilted isobars cause baroclinic pressure gradients and similarly for the pressure gradient in y. Light Coastal systems and open ocean systems absorb light differently for several reasons, so the K P P light and heat fluxes at the surface were altered for the Strait of Georgia model. The parameteriza- tions of light both above and below the surface of the water contribute to the heat and light in the ents. This shows the scheme in ID, the same theory applies for 2D. (2.25) Chapter 2. Methods 39 Strait-. Atmospheric conditions such as cloud cover, time of day and time of year contribute to the solar heat and light budget; attenuation of light once it enters the water is controlled by the clarity of the water. The Strait of Georgia model uses a cloud model alongside theoretical/predicted values for solar parameters to calculate the amount of light that penetrates the atmosphere at a given time and enters the water of the Strait. Several factors control how the light penetrates into the water column, and a final value for irradiance over depth is calculated. Two irradiance profiles result from the light parameterizations: the biologically active part (IPAR), which is used directly in the biological model, and the part that contributes to the heat (JTOTAL), which appears in the non-turbulent heat flux parameterization, QN and contributes to the temperature profile (2.10). There is no biological component in K P P so it was necessary to include an IPAR parameterization either following the examples of the literature, or of our own design. As described below, the traditional incorporation of IPAR w a s deemed inappropriate for the Strait of Georgia, so the latter option was chosen. The inclusion of ITOTAL by Large et al. (1994) is used as a basis for the method followed here, but again, a locally adapted parameterization was developed to best model the Strait of Georgia. Due to the direct impact of the biological part of the light on phytoplankton growth, the scheme to incorporate IPAR is very detailed. However, as the total light is only one of several processes contributing to the modelled temperature, a more generalized process, similar to the Jerlov system used by Large et al. (1994), is employed to calculate ITOTAL- Al l calculations for light use units of W m ~ 2 . To compare with field data (CTD) measured in / ^ E m _ 2 s _ 1 , these units are converted to W m ~ 2 using a conversion factor of 0.2174, recommended by the PAR sensor manufacturer, Biospherical Instruments Inc. (2003). This is based on a quanta to Watt ratio of 2.77x1018quanta s _ 1 per Watt. As the PAR sensor was not calibrated consistently, Chapter 2. Methods 40 this estimate is indeed a rough one, and so surface modelled IPAR cannot be compared to the measured value. The attenuation of light is comparable, however, because this comparison is based on the shape of the irradiance profile rather than absolute values. Light above the water The first parameter we will consider is albedo. Once light from the sun hits the water a certain fraction, the albedo, is reflected back into the atmosphere. A typical albedo value for the earth's surface is 30%; the value used in K P P for Ocean Station Papa is 6% (Large et al., 1994). A value of 18% was calculated for the Strait of Georgia by measuring IPAR just above and just below the surface at several cruise stations on November 11, 2004 (Figure 2.6). A hemisphere plate was attached to the IPAR sensor used in this experiment to prevent re-radiated light from the surface interfering with the measurement of the downward irradiance being measured. The ratio of the measured light below the surface to measured light above the surface was calculated for all the data and the average (18%) was taken to be the new albedo. An "immersion coefficient" of 0.5944 was applied to measurements recorded in air following in- strument calibration to correct for the difference in index of refraction between air and water (Biospherical Instruments Inc., 2003). The data is representative of varied meteorological condi- tions and indicates that 82% of the light measured above the surface remained in the water column at the surface. True albedo varies with season, sea roughness, and water clarity, but is included here as a constant following Large et al. (1994). The cloud model of Dobson and Smith (1988) was used in the K P P model to parameterize the filtration of sunlight through the clouds. This model was deemed less than ideal for the Strait of Georgia because it was developed for the open seas of the North Atlantic Ocean where atmospheric dynamics are significantly different than in coastal regions, so the methodology of Chapter 2. Methods 41 Cloud model Ocean surface Heat flux depends on Qpand P I Jerlov Type 1, 3, or 5 Albedo 18% reflected 44% • 1 PAR Depends on Q F , P, and h Figure 2.5: The light system in the Strait of Georgia. Light begins at the top of the atmosphere with irradiance Qo=1368Wm~2 and is filtered through a cloud model. 18% of the light is reflected, and 82% makes it into the water column. 44% of this light is in the PAR range, and is attenuated depending on river flow, phytoplankton and mixed layer depth; the total light is attenuated according to the Jerlov classification system, which gives a different distribution of light for different values of river flow and phytoplankton. Chapter 2. Methods 42 Figure 2.6: Cartoon depicting experimental set-up for albedo calculation. The C T D equipped with a PAR sensor (the white circle on the top) was held just above the sea surface (1) to measure PAR in the air, and then just below the sea surface (2) to measure how much light was in the water; this was done at each of the nine sampling stations. The ratio of amount of light in the water to amount of light in the air gave an estimate of the albedo (18%), as it told how much light was reflected by the sea surface. Chapter 2. Methods 43 Dobson and Smith was followed to develop a new cloud model for the Strait of Georgia. Dobson and Smith define an atmospheric transmission factor, T as a ratio of downward short- wave radiation, Q, to the mean solar flux, QQ and the sine of the solar elevation, 9 Qo s i n 9 The short-wave radiation data used here is obtained from The Department of Agroecology Science at UBC (Ketler, 2005). Qo is taken to be 1368 W m - 2 , and 9, which represents radiation incident on a horizontal surface above the atmosphere, is a function of time and of the position of the sun relative to the Earth. Cloud fraction data from the Vancouver International Airport (YVR) was categorized according to the okta system that defines a cloud fraction of 9 as an obscured sky and 0 as a clear sky. T can also be defined by the cloud cover T = Ai + Bi sin 9 (2.27) where A and B are coefficients that define how the light is filtered through the clouds, for each cloud type (i). A regression was done on (2.27) and coefficients A and B were calculated for each of the ten cloud categories. The new coefficients (Table A) were smoothed and then replaced the existing K P P coefficients (Table A.2); improvements in model results were noticeable. Figure A.3 shows the improvement of the new fit compared to the existing fit of Dobson and Smith. Surface temperatures were over-predicted using the cloud parameterization of Dobson and Smith (1988) but the new coefficients produced temperatures that matched observations much more closely. Comparison of model calculated IPAR to cloud fraction data (Figure 2.7) shows that the cloud model is working well in a qualitative sense, in that it produces higher light levels when cloud cover is low, and lower light levels when cloud cover is high. The value of incoming radiation at the surface of the water, ITOTAL is defined using these two Chapter 2. Methods 44 03/Feb/02 06 12 15 01/Feb/03 04 07 10 13 16 22 25 28 02/Feb/04 05 11 14 20 23 26 29 02 05 14 17 20 23 26 01/Mar/05 Figure 2.7: Model calculated PAR (line) and cloud fraction (dots) for February of each year (upper panel- 2002, . . . lower panel, 2005) showing that the model does a good job predicting PAR based on cloud fraction: high cloud fraction means obscured sky, and PAR is low. Chapter 2. Methods 45 processes: albedo and cloud filtering. First we replace Qosinff in (2.26) with QA, which represents Qo multiplied by the eccentricity of the Earth's orbit. Substituting (2.27) into (2.26), solving for Q and applying the albedo factor, we get ITOTAL = Qa{A + Bi sin0)(l - albedo). (2.28) The model reads in cloud fraction data (a value, i between 0 and 9) and the appropriate coefficients Ai and Bi are chosen (A.2), giving a value for ITOTAL- Light below the water Light covers a range of wavelengths of 300-2500 nm. Only the light between 400-700 nm can be used for photosynthesis, and is thus termed Photosynthetically Active Radiation (IPAR)- Irradiance of wavelengths greater than 700 nm is termed the infrared (IR). The model calculates the whole light spectrum above the water based on clouds, sun position, and time of day and year. The model requires information about the IPAR part of the spectrum separately from the total spectrum since the IPAR is important for the biological model, and IR is absorbed at the surface as heat which contributes to the parameterization of temperature in the model. Successful parameterization of light at the air-sea interface is necessary to properly model the Strait of Georgia. The K P P model does not differentiate between total light and biologically available light, and employs a method to distribute the light such that all the long wavelength red light is absorbed at the surface, whereas the short blue wavelength light penetrates deeper. This preferential treatment of depth dependent light absorption was replaced in the model with two separate purely exponential decay systems- one for the total light (ITOTAL) and one for the IPAR- The development of each of these schemes are described below. To quantify how IPAR is distributed with depth, a distinction between IPAR and the total light Chapter 2. Methods 46 is required. Integration of the whole light spectrum at 0 m and between 400 and 700 nm indicates that of 44% of total light falls in the 400-700 nm range for clear water (Jerlov, 1976). Surface irradiance, IPAR{® rn) is calculated in the Strait of Georgia model as 44% of the light that has arrived at the surface, ITOTAL(0 m)- Irradiance decays exponentially in a uniform water column at a given wavelength as described by Beer's Law (2.29). To get a depth irradiance profile the attenuation coefficient, K and the surface irradiance value, IPAR(0 ni) are needed. We have already calculated IpAR(Om) above, so it remains only to calculate K. Attenuation is a measure of how light is absorbed in the water and depends on the amount of particulate matter in the water, including organic material such as phytoplankton and inorganic material such as sediment transported by the river. Since the model calculates phytoplankton quantity in the water (P), river flow (Q), and mixed layer depth (h) at every time step, K can be calculated "on the fly" based on these parameters. Equation (2.30) was generated based on P, Q, and h, and calculates a value of K that is used in (2.29) to calculate a depth profile of IPAR, where a\ = 0.011 m 2 m g _ 1 , a2 = 3.92x 1 0 - 5 s m - 4 , 0 : 3 — - 0.012m - 2, and 0 : 4 = 0.27m _ 1. The fit, which has R 2 = 0.41, dictates that as particulate matter in the water column (river sediment or phytoplankton) increases, K increases and light is absorbed more quickly because there is more matter to absorb it. A deep mixed layer depth (/i) indicates that particles are distributed over a greater depth, and light can penetrate deeper into the water column, resulting in a small value for K. Other possible parameterizations for K were investigated but disregarded, including fitting K to P and Q/h (Figure 2.8, R2=0.38), and a fit including only Q and h (R2=0.13); the fit based on Q/h is approximately as good as (2.30) (Figure 2.8c). IpAR(d) = IPAR(0m)e-Kd (2.29) Chapter 2. Methods 47 Apr Jul Oct 2003 Apr Jul Oct 04 Apr Jul Apr Jul Oct 2003 Apr Jul Oct 04 Apr Jul 0 0.1 0.2 0.3 0.4 |Mean K - Fit A| Figure 2.8: Fit choices for irradiance parameterization, (a) Fits A and B fit compared to Mean K, Fit A= fn(P,Q/h), and Fit B=fn(P,Q,/i) (2.30); (b) absolute differences between each fit and Mean K (residuals); (c) Fit A residuals vs Fit B residuals. Fit A and Fit B are approximately equally good fits of Mean K. K = a iP(0 m) + a2Q + a3h + a4 (2.30) The final parameterization for light in the model is for how the total light (300 - 2,500 nm) is absorbed with depth in our system. In the ocean a large part of the incident light is converted to heat. To represent this the model uses ITOTAL to calculate the non-turbulent heat flux profile. The distribution of total light in the water is taken from Jerlov (1976). Jerlov created a system that describes how light is attenuated in the water column by classifying the ocean into several different water types. The main division in his classification system is between open ocean water and coastal water, and within each of these there are smaller divisions based on water properties. Jerlov gives data for each type on how total light is attenuated with depth in the form of a percentage of Chapter 2. Methods 48 total irradiance from sun and sky remaining in the water column at a given depth (1976, Table XXVIII). Only 6 depths ranging from 0 to 20 m are given, and the model needs information on the total light at every half meter from 0 to 40 m, so the Jerlov data was interpolated to fit an exponential distribution with data at every grid point. These data sets are read in by the model at initialization, the model picks a Water Type based on a scheme outlined below, and the total irradiance profile is assigned. The method for determining what Water Types exist in the Strait of Georgia is as follows. Ideally, one could plot total light distribution data for all five Jerlov Water Types and compare it to real data from the Strait of Georgia, and the closest match would determine the Water Type. However, only IPAR data is available from the field sampling program, so a way to classify the water based on these data was required.. Jerlov (1976) gives a 10 m depth averaged attenuation coefficient, K, for several wavelengths in the IPAR spectrum (400-700 nm). A K for each Water Type based only on IPAR and a value for surface irradiance would provide a suitable light profile to compare to the CTD profiles so that Strait of Georgia water could be classified into Jerlov Water Types. To accomplish this, the values for K between 400 and 700 nm were averaged for each Water Type, and this was done based on the distribution of light within that range. It was determined from integration of the whole light spectrum that 45% of I PAR is blue (400-500 nm), 35% is green (500-600 nm), and 20% is red (600-700 nm). Using this distribution a mean K was calculated for each Jerlov Water Type. To compare this to C T D data a 10 m depth averaged K for each CTD cast was calculated. Al l 5 Jerlov Water Type profiles were then plotted using the mean K and compared to each cruise I PAR profile. Each CTD cast was assigned a Water Type based on the Jerlov Type it matched most closely. This analysis of Jerlov Coastal Water Type profiles and S T R A T O G E M irradiance profiles Chapter 2. Methods 49 revealed that Coastal Water Types 1,3, and 5 typically occur in the Strait of Georgia. A seasonal cycle related to the Fraser River runoff and phytoplankton emerged from these comparisons such that during times of low flow and low phytoplankton concentrations (winter) profiles matched Jerlov Type 1. During the peak of the Fraser River flow the water matched Type 3 most often, and at times when phytoplankton concentration is high and Fraser River flow is low the Water Type was most frequently similar to Type 5. To calculate total light distribution with depth, the model needs ITOTAL(® m) and one of the three chosen Jerlov types. A scheme to choose one type based on the seasonal cycle described above was concocted: During times of low flow (Q <3000m 3s - 1) and low phytoplankton concentration (P(0 m)<7mgm~3) the Strait of Georgia water is Jerlov Water Type 1; during the peak of the Fraser River flow (Q >3000m 3s _ 1) the water is Type 3, and when phytoplankton concentration is high (P(0 m)>10mgm~3) and Fraser River flow is low (Q <3000m 3s _ 1) the Water Type is Type 5. The model chooses the Jerlov Water Type and thus the distribution of total light with depth based on Q and P. This scheme picks the correct type 76% of the time, and picks the next best fit 13% of the time. A "correct" choice was defined as times when the fit picked the same Water Type as was chosen by the visual comparison outlined above in which all the data was assigned a Water Type. W i n d Wind is incorporated into the model in the same way as in K P P . Wind speed and direction are measured at 10 m height at Sand Heads Meteorological Station (Figure 1.2) (YVR, 2005) and are transformed from wind speed and direction into meridional (17) and zonal (V) components. Wind stress is calculated according to (A.l) and (A.2) so that the kinematic surface flux of UWQ (2.12) Chapter 2. Methods 50 and VWQ (2.13) can be calculated to be used in the momentum equations (2.8) and (2.9). 2.2.3 B o t t o m Boundar ies The model forces the bottom (40 m) boundary values of T and 5 to be equal to corresponding CTD data. Consider the example of a model run that begins on September 1, 2003, and runs until October 30, 2003, with observations available on October 1. When the model runs to October 1, it reads in 5(40 m) from the CTD file and sets the modelled salinity equal to that value. If there was also a cruise on October 30, the model reads in 5(40 m) from that cruise and the modelled salinity value is updated. This mechanism forces the modelled parameters to be accurate at depth. To do this, at every time step the model searches through the CTD files for one that occurred on the date corresponding to the model's current time stamp. When it finds a match of day and year, it replaces the modelled salinity at 40 m with the observed value. This methodology is admittedly rough, and enough variation in values at 40 m - particularly of N- exists to cause a noticeable change in the modelled water column (Figure 3.14). An inter- polation of the 40 m values between observations would likely produce results that are smoother and more accurate. This modification will be addressed in future studies. When the model is running from 2001 to 2002 there are no C T D data to choose from (see Results chapter for more details on model run specifics), so the code is modified so that the model searches only for a CTD day of year that matches the model date, instead of a model day and year. So for October 1 2001 the model will choose the October 1 2003 C T D file and force the modelled salinity at 40 m to equal the observed value on October 1 2003. Chapter 2. Methods 51 2.3 The Biological Model 2.3.1 Theory The biological model is adapted from the Ecosystem Model of Jeffery (2002). Jeffery's model is an ecosystem model coupled with the K P P model of Large et al. (1994) that models iron-limited phytoplankton at Ocean Station Papa. The model includes zooplankton, a diverse phytoplankton community, and many different sources of nutrients. For use in this project the biological com- ponent of Jeffery's model was simplified to include only one size class of phytoplankton and one source of nutrients to model the particular ecosystem of the Strait of Georgia. The objective of the Strait of Georgia biological model, when coupled with the physical model, is to predict the timing of the spring phytoplankton bloom. As the impact of zooplankton on the arrival time of the bloom is beyond the scope of this thesis, the zooplankton component of the model was replaced with a constant background grazing rate. The size class of phytoplankton used in the model is microplankton, reflecting the species composition in the Strait. The size fraction of the phytoplankton community of the Strait of Georgia greater than 20^m is dominated histori- cally, and was dominated during the STRATOGEM sampling period by diatoms (Harrison et al., 1983; Cassis, 2005). The most common diatom species is the chain forming diatom, Skeletonema costatum which typically dominates the spring bloom (Cassis, 2005; Harrison et al., 1983).. Of sec- ondary importance is the diatom Thalassiosira spp. which appears later in the bloom. S. costatum is a well studied microplankton (Eppley, 1972; Cloern, 1978) as it is a very common species in many of the world's coastal zones. It typically ranges in diameter between 2^m and 20^m and forms chains measuring more than than 20^m. The nutrient most readily available for phyto- plankton uptake in the Strait of Georgia is nitrate (NO 3 ); measured ammonia ( N H 4 ) and urea Chapter 2. Methods 52 Jan Jul 2003 Jul 04 Jul 05 Jul Figure 2.9: Depth averaged (over 0, 5 m) concentrations of chlorophyll (>20^m), silicate, phos- phate, and nitrate (all in pM) measured at S3 during the sampling period. Al l three nutrients are high before the phytoplankton bloom and decline corresponding to the peak of the bloom. concentrations in the Strait are trace (Harrison et al., 1983). Silicate and phosphate concentration at S3 is high before the bloom and decrease with the peak of the bloom (Figure 2.9), similarly to N. Phytoplankton require only one mole of phosphate and 15 moles of silicate for every 16 moles of nitrate required (based on the Redfield ratio of N:Si:P = 16:15:1). Examination of nutrient concentrations preceeding the blooms of 2003, 2004 and 2005 show that each year nitrate will become the first limiting nutrient once the bloom begins ([A7]/16 < [P]/l < [Si]/15). Therefore nitrate was the only nutrient included in the model. The model was adapted in such a way that the zooplankton component, the other size classes and the other nutrient sources can easily be added back in for further studies. Chapter 2. Methods 53 Biological Model Numerics The growth rate of phytoplankton and rate of change of nutrient concentration are calculated based on Jeffery (2002) with the following terms excluded: an active zooplankton component, preferential treatment of nitrate, detritus and particulate organic nitrogen components. Table B . l gives values for all biological parameters. The biological model is a Nitrogen model, that is, all the biological units are nitrogen. The phytoplankton in the model is initialized using CTD fluorescence measurements, which are con- verted from mg m - 3 chl to LtM Nitrogen using a conversion factor of 1:1. This ratio was arrived at by assuming a gCarbon : gchl ratio of 80:1 (Parsons et al., 1984) and a Carbon : Nitrogen molar ratio of 6.625:1 (Redfield). Nutrients are initialized using sampled nitrate (N) bottle data. The data from the bottles (0, 5 , 10, and 30 m) is interpolated at half-meter intervals to match the grid spacing of the physical model. N concentrations between 30 m and 40 m are assumed to be uniform; this is consistent with observations at other stations in the Strait for which deeper bottle samples are taken. When nutrients are not limiting phytoplankton growth, uc, depends on light availability " - ' ^ • M ^ ^ X r r / ^ ) <2-31> where P growth increases with light with an initial slope of o until saturation, or maximum growth {Rmax), is reached. If light continues to increase above levels of saturation, growth inhibition occurs ((3). Growth can be slowed by losses due to exudation ( 7 ) , or excretory carbon losses. When light is not limiting phytoplankton growth depends on nutrient availability R N i<-maxly i n n < - > \ un = (2.32) where K is the half-saturation uptake constant for nitrate (./V). Phytoplankton in the Strait of Chapter 2. Methods 54 Georgia take up nitrate during dark periods (night) (Price et al., 1985). However, this process cannot be modelled here because dark conditions imply light limitation, and therefore no nutrient uptake can occur. A temperature dependency of maximum growth rate is implemented in the model following Eppley (1972) using R 0 R-max = Ro * QlO ^ * (2.33) where Tref is the optimum growing temperature (20°C), Qio=l-88 (Eppley, 1972), and T is mod- elled temperature. This function is executed at all model depths. The maximum growth value for S. costatum at optimum growing temperature, R 0 , and values for a and (5 in (2.31) were determined experimentally by fitting them to the growth data for S. costatum of Jorgensen (1966) represented as a growth versus light curve in Cloern (1978). Phytoplankton growth (2.34) and nutrient consumption (2.35) depend on the availability of IPAR and N and are calculated as a function of whichever is minimum. Growth is also limited by losses due to grazing and non-temperature dependent phytoplankton respiration, R m . Diffusion of P and N is the same as the diffusivity for salinity, Ks(a), (2.3), and upwelling W follows (A.20). ^ = jminMW), un(N)} - R j j P + A Jtf a|P J + w { p ) ( 2 3 4 ) ~ = -!(mm[uC(IPAR),un(N)}}p+ A j ^ ^ j + W(N) (2.35) Mortality (zooplankton grazing) is included in the general loss term, R m , a linear function of phytoplankton biomass, (P x R m ) , rather than a specific grazing function because it is an unknown parameter to which the model is very sensitive, as will be discussed in the M o d e l Tuning and Sensitivities of the Results chapter. A grazing function would therefore only Chapter 2. Methods 55 add to the uncertainty surrounding the zooplankton grazing. Rm was determined by trial and error using values from the literature as a starting point. The value chosen for Rm (Table B. l ) produced optimal model results; the Mode l Tuning and Sensitivities section defines "optimal" and discusses this choice. Equations (2.34) and (2.35) are solved using an implicit-explicit scheme (B.2) wherein the first terms in (2.34) and (2.35) are solved using a fifth-order Runge-Kutta with adaptive timestep, as detailed in Appendix B . 2.3.2 Bottom Boundaries The bottom boundary condition used for the physical parameters outlined in the Physical M e t h - ods section also applies to bottom boundary conditions for P and N. This ensures that conditions at 40 m depth always match observations; this is especially important for nitrate, as the forcing ensures that levels of at the bottom stay high. If N was allowed to deplete at depth, the water that is upwelled to the surface would be A^-depleted and P would have low food supply. 5 6 Chapter 3 Results 3.1 Physical Model Results The physical model was run independently of the biological model over various time scales and compared to field data. Model parameters were analyzed over depth and time and adjusted so that results were representative of the Strait of Georgia. In addition to C T D data, which was used to verify model output in depth, ferry data was used to verify parameters at the surface over time. The physical model was considered to be working properly when model profile output matched the general shape of CTD profiles. In particular, the mixed layer depth, surface property values, and depth and steepness of pycnocline, halocline or thermocline must match. As the values at 40 m were forced to be equal to CTD values at 40 m (see Methods chapter); profiles below 30 m are generally well matched. 3.1.1 M o d e l T u n i n g and Sensi t ivi t ies The physical model is sensitive to the internal mixing parameters, specifically the parameters con- trolling interior shear mixing (VQ) and unresolved internal wave breaking (v™). In KPP , Large et al. (1994) choose VQ based on observed background diffusivities, but data for the Strait of Georgia is largely unknown so the value is chosen by trial and error. Increasing the diffusivity results in Chapter 3. Results 57 increased mixing and a deeper mixed layer compared to a smaller value of diffusivity which causes less mixing and a shallower mixed layer. The same trend is true for the internal wave shear values (both a scalar and a momentum quantity, and v™m), so that a higher value increases mixing and deepens the mixed layer depth, and so forth. The chosen values result in well-matched profiles under'varying conditions throughout the year; ultimately only the value for the maximum internal shear VQ is changed from the original (from 5 x l 0 ~ 3 m 2 s _ 1 to 1 x 10~ 3 m 2 s _ 1 ) . Upwelling (vertical advection) is considered the most important parameter in terms of matching modelled salinity and temperature profiles to observations, as the model is very sensitive to changes in the value. The slope of the pycnocline, the mixed layer depth, and surface values can all change drastically if the upwelling value is changed. Specifically, increased upwelling causes surface salinity to increase and surface temperature to decrease with the addition of cold, saline water to the surface from below; surface salinity decreases and surface temperature increases with weakened upwelling. Stronger upwelling causes the mixed layer to shallow and weakened upwelling causes the mixed layer to deepen. The area over which the upwelling occurs is the fraction of the concentric circles, 1 /Cw (Figure A.2). The optimal value for this variable will result in an ideal upwelling situation, where surface values of physical properties are predicted correctly by the model. The model is quite sensitive to different values of Cyv, variations of which represent different upwelling strengths. It seemed that a seasonal cycle of upwelling correlated strongly to the Fraser River and related to the turbulence of the flow (Figure 1.8), and it was thought that implementing this seasonal pattern into the upwelling scheme instead of using one value at all times would improve model performance. It was found that the turbulence reaches a maximum typically a few days before the flow begins to increase in March or April; at this time we see the greatest vertical transport of water from Chapter 3. Results 58 the lower layer into the upper layer (Riche, 2005). During times of low flow and corresponding low turbulence, the vertical velocity of transport is on a scale of 5 x l 0 5 m s _ 1 less, or half as fast as during the peak turbulence. These values are significant and support the idea of a seasonal cycle scheme for upwelling in the model. However, testing various values of Cw in the model and various seasonal schemes did not ultimately provide an optimal upwelling solution; rather, the optimal solution was found to be CV=3, indicating that upwelling is taking place over an area one third the size of the area defined by the circles in Figure A. 2. 3.1.2 Vertical Variability The figures presented in this section and the Biological Vert ical Variabi l i ty section were chosen to represent a variety of seasonal conditions, and a balance of good and poor results. For each of the figures (Figures 3.1- 3.5) the model was initialized with CTD data from a given cruise and run until a later cruise so that the model results could be compared to observations. A run initialized with the May 9, 2002 cruise C T D data, and run until the next cruise on June 6, 2002,(Figure 3.1), and a second run from April 23 to May 28, 2003 (Figure 3.3) show different results under similar high Fraser River freshet conditions. Winter conditions are illustrated by a run initialized with Dec.4, 2002 cruise C T D data, and run until the next cruise on January 10, 2003 (Figure 3.2); whereas late summer conditions are illustrated by a run from July 3 to Aug.13, 2003 (Figure 3.4). The run initialized on March 2 to March 17, 2004 represents results both from a shorter model run and a run during spring bloom conditions (Figure 3.5) because during the spring bloom cruises occurred at a higher frequency. Chapter 3. Results 59 Temperature The model predicts temperature consistently well. Figures 3.5 and 3.2 show modelled surface tem- perature is within 0.5°C of observed surface temperature, and Figures 3.3 and 3.4 show modelled and observed surface temperatures that are within less than 2.4°C of each other. Figure 3.3 is an example of the well-matched modelled and observed thermocline and well-matched profiles below the mixed layer depth. Salinity Salinity is predicted well by the model, but less consistently than temperature. Figures 3.2 (on a different scale than the others) and 3.4 show well-matched modelled and observed surface salinity, with surface values within 1 and 2 of each other, respectively. Figures 3.1 and 3.3 are examples of poorly-matched surface salinity, with modelled and observed values at the surface differing by more than 6 and 8, respectively. These both represent high Fraser River flow conditions, and are perhaps an indication that the model does not do as well at predicting very low salinity. The halocline is modelled well in most cases (Figures 3.2, 3.3) although there are also cases in which the modelled halocline is too shallow (Figure 3.4) or too deep. The small scale (zoom) of Figure 3.2 shows a well-matched halocline and well-matched variations below the mixed layer. Density The density profiles mirror the salinity profiles very closely in most cases (Figures 3.2 through 3.5). Figure 3.1 shows that although the surface values of the modelled and C T D estimated density are not well-matched, the slope of the pycnocline is very well-matched. Chapter 3. Results 60 Irradiance Irradiance profiles are predicted moderately well by the model, which computes the attenuation coefficient (K) based on Q> P, and h (2.30). Figures 3.2, 3.3, and 3.5 are examples of a well predicted K, as the slope of the modelled light curves are matched very well to the observed curves; other examples show poor prediction of attenuation slope (Figures 3.1 and 3.4). Absolute values of modelled and measured surface irradiance are incomparable due to the calibration of the PAR sensor, but it is valid to compare the attenuation of the two curves. Eddy Viscosity and Diffusivity Eddy viscosity and eddy diffusivity are plotted to show the relationship between the physical pro- cesses governing the mixed layer and the lower layer. The mixed layer depth is defined in the physical model as the depth to which eddies will penetrate, represented visually by the relation- ship between eddy viscosity, eddy diffusivity, and density (Figure 3.6). Eddy diffusivity (scalar properties: temperature T, salinity S) and eddy viscosity (velocities: u, v) (2.3) depend on the mixed layer depth, h, the turbulent vertical velocity scale, wx, and the the vertical shape function, G(cr). Al l these values are defined so that they are the same for temperature and salinity (see Methods chapter), so eddy diffusivity profiles for T and S are identical and are plotted as Eddy Diffusivity (Figure 3.6). Eddy viscosity and diffusivity are nearly identical in profile shape and magnitude, with diffusivity greater than viscosity in most cases. Eddy viscosity and diffusivity are maximum in the upper mixed layer and quickly reduces to almost zero at the mixed layer depth (Figure 3.6 (a) and (c)). Small increases in eddy viscosity below 35 m (Figure 3.6 (a) and (c)) coincide with a slight change in stratification not visible in the density profile at this scale. Stronger mixing in panels (a) and (c) result in a deeper mixed Chapter 3. Results 61 Salinity 25 - >T 30 - , J 35 - \ 401 l_j 1 5 10 15 20 25 Temperature f C) Figure 3.1: (a) salinity (dashed lines), temperature (solid lines); (b) light (dashed lines), density (solid lines). Model predicted (thick lines) and observed data (thin lines). May 9- June 6, 2002. Chapter 3. Results 62 Figure 3.2: (a) salinity (dashed lines), temperature (solid lines); (b) light (dashed lines), density (solid lines). Model predicted (thick lines) and observed data (thin lines). December 4 2002- January 10 2003. NOTE: Temperature and salinity are plotted on a different scale than in the other figures. Chapter 3. Results 63 Figure 3.3: (a) salinity (dashed lines), temperature (solid lines); (b) light (dashed lines), density (solid lines). Model predicted (thick lines) and observed data (thin lines). April 23- May 28, 2003. Chapter 3. Results 64 PAR (Wrrf •=) 100 200 10 15 20 25 Temperature f C) Density (kgm Figure 3.4: (a) salinity (dashed lines), temperature (solid lines); (b) light (dashed lines), density (solid lines). Model predicted (thick lines) and observed data (thin lines). July 3- August 13, 2003. Chapter 3. Results 65 Salinity PAR (Wrrf 2) Temperature i°C) Density (kgm 3 ) Figure 3.5: (a) salinity (dashed lines), temperature (solid lines); (b) light (dashed lines), density (solid lines). Model predicted (thick lines) and observed data (thin lines). March 2- March 17, 2004. Chapter 3. Results 66 Density (kgm" 3 ) Density (kgm 3 ) Density (kgm" 3 ) Figure 3.6: Modelled eddy diffusivity (m2s : ) , eddy viscosity (m2s and density (kgm 3 ) . (a) March 8 to April 7 2005; (b) May 8 to June 5, 2002 (b); (c) March 7 to March 19, 2003. layer than in panel (b) because stronger eddies penetrate deeper. In panel (b) diffusivity is low, eddies are weaker, and so the mixed layer is shallow. 3.1.3 Temporal Variability Valuable information and understanding of how well the model is working can be attained from plotting model output throughout the entire model run. Analysis of the response of modelled parameters to different physical forcings can be discerned and insight into the performance of the model can be gained. In this section several plots are presented that show the interaction between Chapter 3. Results 67 physical parameters modelled over long time periods. These results show that mixed layer depth is controlled by wind until the arrival of the Fraser River freshet, which stratifies the water column and decreases the mixed layer depth. This stratification is sometimes disrupted by wind events. Wind Wind is the sole mechanism through which kinetic energy is introduced into the model and its effect on the mixed layer depth is clear (Figure 3.7). The mixed layer depth is highly dependent on the wind, as we see the deepest mixed layer depths coinciding with strong winds, especially in the winter of 2002. A wind event in late August 2003 deepens the stable, shallow mixed layer from less than 5 m to greater than 30 m; the same type of event occurs in late May 2005. Fraser River The biggest source for potential energy input into the Strait of Georgia is the Fraser River. The effect of the Fraser River freshet on the mixed layer depth is to stratify the water column and decrease the mixed layer depth (Figure 3.7). The mixed layer depth is unstable and frequently deeper than 20 m in the fall and winter, but shallows to less than 10 m when the Fraser River freshet arrives. Wind events can disrupt the stratification produced by the layer of fresh water (late August 2003) but the freshwater flux is the most significant factor determining the mixed layer depth during times of high Fraser River flow. Surface salinity is clearly controlled by the Fraser River (Figure 3.8), an indication of a suc- cessful freshwater flux parameterization (2.21). The effect of the Fraser River is mainly limited to the top 10 m of the water column, with the strongest freshwater signal confined to the top 5 m (Figure 3.9). Chapter 3. Results 68 Figure 3.7: Modelled mixed layer depth, wind magnitude cubed, and Fraser River flow for Novem- ber 1 to October 1 for each year. Strong winds deepen the mixed layer depth, while weak winds are associated with shallow mixed layer depths. The freshet coincides with a stable, shallow mixed layer depth. Chapter 3. Results 69 Temperature The resolution of diurnal temperature cycling (Figure 3.10) indicates that the model can resolve small-scale effects as well as large-scale effects such as those produced by the Fraser River. Temper- ature fluctuations in the Strait of Georgia are caused by air temperature and heat from sunlight, and a strong correlation exists between air temperature and surface temperature; even a dramatic shift in air temperature (Figure 3.10b) is well-modelled. Wind mixing can cause deeper, colder water to be mixed to the surface (Figure 3.10a). The annual cycle of temperature is also evident in model results (Figure 3.9). Water up to 40 m depth is warmed during the summer and cooled during the winter, but the coolest winter water occurs at the surface. The warm summer surface water can extend down to 10 m. Clouds and Light Model calculated light compares well to buoy data (Figure 3.11) collected at the Halibut Bank Buoy (Figure 1.2) (Gower, 2005). The buoy measures IPAR at the surface above the water (in the air), so the albedo is applied to the data to make it comparable to modelled surface IPAR, which is calculated for the surface of the water (underwater). The good correlation between modelled and observed surface irradiance (Figure 3.11a) is evidence that the estimation of light at 0 m done by the model, particularly the estimation of how much of the total light is IPAR, is correct. Comparison to the buoy data under cloudier conditions (Figure 3.11b) is a better test of the cloud model, and we see that although there is a slight mismatch, overall the cloud model parameterizes clouds well. This mismatch is undoubtedly due to the distance between Halibut Buoy and Y V R , the location of the cloud monitoring station (Figure 1.2), resulting in different cloud cover at different times. Chapter 3. Results 70 Figure 3.8: Model computed surface salinity (thick line) and Fraser River flow (thin line). Modelled salinity responds well to the Fraser River, decreasing with increased flow and decreasing when flow is high. Chapter 3. Results 71 Figure 3.9: Contours of modelled salinity (upper panel) and temperature (lower panel) from September 2001 to July 2005. The freshwater signal of the Fraser River freshet appears in the spring of each year; the seasonal cycle of temperature is also visible. Chapter 3. Results 72 30 9_ 25 Q. I 2 0 3 15 10 » 2000 E a 1000 I I I 1 I I I I i i i i 1 1 1 — Air temperature | * — - Modelled surface temperature i l i i i i Jb) r \ \ , A . A . . , J * \ . . r v \ . f \ As £ IO/Jul/03 12 to E 2000r •D I 1000 0 I I J") 1 i Ik i i . . . . . - . - A 1 & L. i .1/ |j. « i « i 15 o o IOM<» ? 5 03/Jan/05 08 13 02/Feb Figure 3.10: Model computed surface temperature and air temperature at Y V R for July 10 -31, 2002 (a) and Jan. 1 - Feb.20, 2005 (d); wind magnitude cubed (b) and (c) during the same time period as (a) and (d). The cooling of the surface temperature on July 22 and 28 is not caused by the air temperature, which does not dramatically decrease, but is caused by wind events on the same days (b) which cause deeper, colder water to be mixed to the surface. Surface temperature variations on January 16-28, 2005 are due mainly to air temperature variations (c), (d), as the curves are closely matched during this period. The wind event on January 13 causes warmer water from below the surface to be mixed to the surface. Chapter 3. Results 73 3501 1 1 1 r Figure 3.11: Modelled IPAR(0 m) (dotted line) compared to observed IPAR(0 m) (solid line) mea- sured at Halibut Bank Buoy (Figure 1.2) for (a) July 27 - August 5, 2003, (b) Oct. 25- Nov.16, 2003. Panel (a) is during a relatively cloud free period and shows that the model does a good job predicting IPAR- Panel (b) shows that during cloudier conditions the model sometimes over-predicts and sometimes under-predicts I PAR, but overall does a good job representing cloud conditions. Chapter 3. Results 74 S e a s o n a l V a r i a b i l i t y Model runs over seasonal time scales, from October to March (winter), March to May (spring), and May or June to September (summer) were analyzed to assess possible seasonal consistencies or differences. The results are as follows. The model does a better job at predicting differences between modelled and observed properties at the surface in the winter than in other seasons, as there is less variability in salinity and temperature. Modelled surface temperatures are within 0.4°C of the data, modelled surface salinity values are within 2, and modelled density surface values are within 2 k g m - 3 of the data. Mixed layer fluctuations in the winter are mostly due to wind (Figure 3.7), the treatment of which is handled well by the physical model, as it predicts mixed layer depth within 2 m for each winter run. Slopes of the pycnocline, halocline and thermocline calculated by the model are in excellent agreement with the data during the winter months. Mixed layer depth in the observational data is estimated as the depth where consecutive densities differ by more than O.Olkgm - 3 . Temperatures are typically over-predicted at the surface for the spring runs, resulting in water that is from 2 to 4.5°C too warm. Thermoclines are generally the same slope and depth as the CTD data, and mixed layer depths are very well predicted. Salinities at the surface are typically more saline than the field data, indicating that the modelled freshwater signal of the Fraser River is too weak, i.e.: not enough freshwater in the model. Halocline slopes are all good, and mixed layer depths are predicted very well. General shape of profile and slope of pycnocline are excellent. Temperatures at the surface during the summer runs match data well, with some runs slightly too cold and some runs slightly too warm. Mixed layer depths are good, but thermoclines are often too flat and shallow. Salinity profiles are often somewhat too saline at the surface, but haloclines have the right slope and mixed layer depth is good. Chapter 3. Results 75 3.1.4 Tes t ing M o d e l Performance The way in which the model responds to varied conditions in meteorological forcing such as high winds or increased freshwater input is a good indication of its success in modelling the physical system of the Strait of Georgia. D a t a F o r c i n g Running the model with winds forced to be 11 m s - 1 (equivalent to 40kmh _ 1) for an entire run, regardless of run length or season, causes a deepening of the mixed layer and pycnocline from the extra input of kinetic energy and extra mixing. Forcing winds to be <1 m s - 1 (2kmh _ 1) causes the mixed layer and pycnocline to shallow due to the absence of mixing. Low Fraser River flow (set to l x l O 3 m 3 s _ 1 ) causes the surface salinity to increase significantly, e.g. by 10. Light is able to penetrate much deeper when the sediment-laden Fraser River flow is low (lower attenuation coefficient). When Fraser River flow is set to a constant of 5 x l 0 3 m 3 s _ 1 for a whole season the surface salinity decreases and light is absorbed at shallower depths (higher attenuation coefficient). The model responds appropriately to forced clear skies as well as forced cloudy skies. Running the model with clear skies without changing the humidity or air temperature for an entire run resulted in warmer temperatures at the surface and at depth, regardless of model run length or season. Light at the surface is consistently higher than with true cloud conditions, and also slightly higher at depth. Running the model with cloudy skies decreases the modelled irradiance, and decreases surface temperatures, as expected. % Chapter 3. Results 76 Ferry Comparisons Comparing model output to ferry data provides great insight into what is happening in the Strait of Georgia over large temporal and spatial scales. The discrete profile provided by a CTD cast lends no clues to the past or present, or to the water adjacent to S3, and bottle sampling can miss important events in the water column, so ferry data can be a valuable comparison tool. The ferry data that corresponds to S3 is calculated by extracting ferry data that falls into a box around S3 defined by one minute of latitude and longitude on either side of S3 (49.06:49.08, -123.34.-123.32). A thermosalinograph on the ferry measures the surface water that is pumped through the ferry system to be used as engine coolant. This water is sampled from ~3 m depth near the back of the ferry, but often the sampled water has undergone mixing caused by turbulence created by the ferry, so the actual depth of sampled water can be anywhere in the top 5 m of the water column (Halverson et al., 2005). Temperature fluctuations measured by the ferry match the modelled surface temperature well (Figure 3.12), with temperatures throughout the winter matching to within 1°C. The model typi- cally over-predicts temperature in the summer by between 1 and 2°C. A comparison of salinity is less straightforward (Figure 3.13) because the ferry salinity record fluctuates on a greater range than ferry temperature record. The fluctuations are thought to be due to the influence of the tide and the wind as it pushes the plume over S3 and then away from S3 several times over the course of a day. Since the ferry records data every second, the salinity record for one day can range from 30 to 10. Despite the disagreement in absolute values of salinity (some of which are due to equipment calibration issues (Halverson et al., 2005), it is clear that the model does a good job depicting the overall salinity. Small-scale salinity events such as heavy rainfall (October, 2003 and January, 2005) are apparent, as well as the larger-scale fluctuations. Figure 3.12: Modelled 3 m temperature (black line) and temperature measured by the ferry at ~3 m (grey dots) for October - November of 2003, 2004, and 2005. Chapter 3. Results 78 Figure 3.13: Modelled 3 m salinity (black line) and salinity measured by the ferry ~3 m (grey dots) for October - November 2003, 2004, and 2005. Modelled salinity is somewhat over-predicted, but large-scale variations are temporally matched. Chapter 3. Results 79 3.2 Biological Model Results The Strait of Georgia coupled biological-physical model was run for the entire sampling period, one year at a time. Except for the first year of sampled data, which began in April 2002, each model run began in October and the model ran for 13 months (Table B.2). The unknown biological model parameter, Rm, was tuned to predict the spring bloom when initial conditions were mid-autumn. A different value is required if the model is initialized at the peak of the bloom, or in the summer. For this reason, the model was not run beginning with the first S T R A T O G E M data in April 2002, as results were inconsistent with runs beginning in October. To predict the 2002 bloom, a CTD profile from September 2001 (Linguanti, 2001) was used to initialize the model. Al l forcing data is available for 2001, but bottom boundary conditions are forced with data from other years, as described in the Methods chapter. The physical model was not adjusted to match biological parameters, as the implementation of the biological model, aside from the dependence of the light parameterization on P, does not affect the running of the physical model. 3.2.1 M o d e l T u n i n g and Sensi t ivi t ies Leading up to the spring bloom, the coupled biophysical model is not sensitive to the half-saturation constant for nitrate, K, because this value is so small (Table B. l ) and nitrate is plentiful- much greater than K. The model is relatively insensitive to 8, the light inhibition constant, and cr, the initial slope of the growth-due-to-light curve. A drastic increase in a (>50%) can cause the bloom to arrive earlier, and an equally unrealistic decrease can delay the bloom. The magnitude of the bloom decreases slightly with a larger value for 3. The model is most sensitive to changes in the value for growth losses, Rm, the only parameter Chapter 3. Results 80 to which the biological model was tuned. In general, increasing losses decreases bloom magnitude and delays the arrival of the bloom, while decreasing losses increases bloom magnitude and causes it to occur earlier. Small variations of the losses by <1% can change the bloom arrival time by as much as one week. The great variability due to the loss term (Rm) caused considerable difficulty in obtaining a set of parameters that could be trusted. An ideal set of parameters would change the bloom arrival time gradually as parameters were changed gradually, but this was not always the case. Often small changes in the loss term caused large bloom arrival delays, while other times large changes in the loss term had no effect on the arrival time. In one scenario the loss term that produced the best results for all four years was unreliable as a small change (0.01 x l O - 6 s _ 1 ) in losses caused the bloom in 2005 to arrive 10 days earlier than with the previous value of Rm. This was an unfavourable and un-physical result, and forcing tests did not perform as expected (e.g. no cloud cover caused a later bloom than total cloud cover). This indicated that the set of parameters was invalid and therefore it was discarded. Although the model was not sensitive to the parameter for nitrate uptake (n), the surface nitrate was affected by changes in the respiration and mortality (the loss term), which often produced large peaks in the modelled surface nitrate that did not correspond to observations. The parameters that were ultimately chosen produce a minimum amount of these features. The loss term that was ultimately chosen produces physical and logical results, but results in modelled surface P that is lower than observations. A lower value of Rm increases P, but also causes the bloom to arrive too early. Since emphasis is on matching bloom time arrival rather than matching P magnitude, the value for Rm that gives the best timing is chosen, not the value that gives the higher P. Different values of Rm produce a more realistic replenishment of surface Chapter 3. Results 81 nitrate in the fall but often at the cost of the bloom time arrival. The bloom arrival time is the most important aspect to model correctly, so the values that achieved that were the values chosen. 3.2.2 Vertical Variability The model produces good results at depth as well as at the surface (Figure 3.14). Chlorophyll concentrations during the bloom extend to depths of more than 10 m below surface at times, and the abundance of nitrate throughout the year at depth is represented well. Nitrate is depleted to the same depths to which the phytoplankton bloomed, and this depletion occurs following a bloom event. Another measure of the success of the model is the comparison of the modelled biological properties to the CTD observations. Each panel of Figure 3.15 shows modelled P and N and observed fluorescence (CTD) and nitrate (bottle). From left to right the panels in Figure 3.15 are results from model runs that span the following dates: (a) April 9 to April 16, 2003; (b) February 23 to March 8, 2005; (c) October 8 to November 5, 2002. Panel (a) represents the peak of the bloom (b) represents beginning of the bloom; (c) represents fall conditions. The model typically under-predicts surface phytoplankton (Figure 3.15 (b), (c)) compared to fluorescence because there is only one phytoplankton group, and one phytoplankton size fraction (>20pm) included in the model. CTD estimated fluorescence measures total chlorophyll in the water column- all size fractions, plus anything else in the water column that might fluoresce. For this reason, surface modelled values often match bottle sampled chlorophyll (>20p,m) better than they match CTD fluorescence (Figure 3.15 (b), (c)). The >20pm size fraction makes up, on average, 50% of the total chlorophyll (>20/um + >2/j,m + >0.2/^m) for all S T R A T O G E M data at S3 at 0 m. For the model runs represented in panels (a) and (b) the >20pm fraction makes Chapter 3. Results 82 up almost 100% of the total chlorophyll. As these represent conditions at some point during the bloom, this domination by large phytoplankton is expected. During the low biomass conditions represented in panel (c), we expect smaller phytoplankton to dominate. Indeed, the >20pm size fraction make up only 25% of total chlorophyll. Comparing modelled phytoplankton on this date to 25% of the fluorescence results in a value closer to 3.75 m g m - 3 , much closer to the modelled value. This implies that, in the case of winter conditions, the mismatch between modelled phytoplankton and fluorescence can be accounted for by the absence of smaller phytoplankton. This assumption is not valid during bloom conditions. The model often correctly predicts the slope of P depletion with depth as well as the depth of the subsurface chlorophyll maximum (Figure 3.15 (a), (b)). Nitrate is well matched in panels (a) and (b) at the surface, but is under-predicted between 0 and 5 m. These low levels of nitrate compared to measured values could also contribute to the low modelled phytoplankton. The slope of N depletion with depth is well predicted (a) and (b); (c) does a good job reproducing the mixed layer depth but greatly under-predicts nitrate at the surface. The observed nitrate data is interpolated bottle data, for which true data points only exist at 0, 5, 10 and 30 m depth (open circles), so the modelled nitrate is often smoother (b) than observed, but may have some more curves that the discrete bottle sampling may have missed (a). 3.2.3 Tempora l V a r i a b i l i t y The arrival of the spring phytoplankton bloom in the Strait of Georgia is predicted successfully by the coupled biophysical model (Figures 3.16 and 3.17). In each of the four sampled years (R l through R4, Table B.2), both the phytoplankton bloom arrival, peak, and decline, and the nitrate depletion occur at the same time as is observed in the field. Nitrate levels throughout the winter increase to levels that are seen in the field, while P decreases to low levels during the light limited Chapter 3. Results 83 Figure 3.14: Contours of modelled log of chlorophyll (upper panel) and nitrate (lower panel), both in mg m _ 3 c / i / , from September 2001 to June 2005. The phytoplankton bloom extends to more than 10 m depth; surface nitrate depletion is extensive following a bloom event. Chapter 3. Results 84 Phytoplankton (mg m chl) 10 20 „ 30 40 Phytoplankton (mg m~ chl) 10 20 30 40 10 20 30 40 Nitrate (uM N) Phytoplankton (mg m chl) 10„ 20 30 40 0 10 20 30 Nitrate (uM N) 10 20 30 Nitrate (uM N) Figure 3.15: Phytoplankton (solid lines) and nitrate (dashed lines) observed (thin lines) and mod- elled (thick lines). Observed P is CTD fluorescence, (a) April 9 to April 16, 2003; (b) February 23 to March 8, 2005; (c) October 8 to November 5, 2002. Open circles are bottle nitrate samples; asterisks are bottle chlorophyll samples greater than 20 nm. Chapter 3. Results 85 winter months. Fall blooms are also observed and modelled accurately in September 2003 and 2004. The second and third peaks of the 2005 bloom are well modelled. The model does not do a good job estimating the magnitude of the blooms observed in the data, but this is not unreasonable because of the simplifications made to the biological model. Flagellates are present in the Strait of Georgia year round, and low modelled P can be explained in part by their absence and in part by the constant loss term, Rm, which is neither temperature dependent nor does it include a seasonal effect for phytoplankton mortality that might represent zooplankton grazing. The fall increase in Nitrate does not occur with the chosen biological parameters but can appear with a different value for Rm. The nitrate is forced at 40 m every cruise to match CTD data (see Methods chapter) which introduces high levels of N to the bottom, but upwelling does not overcome the stratification that exists in the modelled results (Figures 3.9, and 3.14). This stratification does not agree with stratification that appears in the observations (Figures 1.3 and 1.7). This is a matter that warrants further investigation but as the spring bloom is the focus of this thesis, the fall nitrate renewal will have to be left for future work. 3.2.4 Tes t ing M o d e l Performance Data Forcing The performance of the biophysical model can be assessed by running the model with artificial data, or manipulating the modelled physical processes to ensure they are working properly. The process of upwelling should distribute nutrients and phytoplankton throughout the water column so that nitrate, which exists in high concentrations at depth, is continually being supplied to the upper layer for biological uptake. Upwelling should also ensure that phytoplankton populations, Chapter 3. Results 86 Figure 3.16: Modelled surface phytoplankton for all four modelled years (solid lines) and greater than 20 pm bottle sampled chlorophyll (asterisks). Chapter 3. Results 87 Figure 3.17: Modelled surface nitrate for all four modelled years (solid lines) and bottle sampled nitrate (dotted lines). Chapter 3. Results 88 which are greater near the surface than at depth, are continually being advected or diffused away. With upwelling set to zero in the model over a short time scale phytoplankton populations at the surface initially increase as expected until surface nitrate runs out and has no mechanism for renewal. Maintaining these conditions over a longer time scale causes P to decrease due to lack of food and extra competition for the little available food. Increased nutrients and salinity in the mixed layer can be achieved in the Strait of Georgia through increased entrainment or increased wind. Running the model with artificially forced high upwelling (w—w x 10) causes an increase in N in the surface layer, providing food for P and causing an increase in their population. Analysis of transport rates of N and P reveals that nutrient flux (in f i M m " 2 ) is two orders of magnitude greater than P flux (in ^ M m - 2 ) because there are more nutrients introduced at the bottom to be upwelled away. Therefore, although there is increased P loss to the surface layer when upwelling is high due to advection the net result of P is to increase due to the increased N. Strong winds deepen the mixed layer and mix phytoplankton through the layer to the bottom of the mixed layer depth (Figure 3.18). This is apparent by the decrease in P in the top'20 m and the increase in P below 20 m (Figure 3.18b). Total P biomass declines for a brief period (bottom panel, Figure 3.18a) because most of the population is now deep in the light limited zone and growth is slowed by this lack of light. If growth slows enough so that loss terms are greater than growth terms (2.31), P populations decrease, and if the population decreases below zero the model sets P=0. The depth at which this occurs is called the critical depth (Figure 3.19). It is clear that P that become mixed below the critical depth decrease in biomass. These results provide strong evidence that the biological model is responding properly to physical forcings. Figure 3.18: Model results from March 2003 showing the effect of strong winds, (a) March 1 to March 15. From top to bottom: wind magnitude cubed; mixed layer depth; modelled chlorophyll at 0,3,6,10,15,20,30, and 40 m depths; depth integrated modelled P. (b) Details for first wind event, March 3. Wind magnitude cubed, mixed layer depth, and modelled P (mg m~ 3 chl) depth contours. Chapter 3. Results 90 Figure 3.19: Mixed layer depth, modelled 0 m phytoplankton, and critical depth for Feb.l to March 27, 2003. When the mixed layer depth is deeper than the critical depth there is a decrease in phytoplankton. Chapter 3. Results 91 Figure 3.20: Modelled 3 m phytoplankton (black lines) and ferry recorded chlorophyll (grey dots) for October to November 2003, 2004, and 2005. Model results for spring bloom arrival agree well with ferry records, within days of each other for each case. Chapter 3. Results 92 Ferry Sampling Ferry measured chlorophyll is one of the best records of the bloom dynamics in the Strait of Georgia, enabling us to fill in the gaps that are created by discrete C T D and bottle samples. The ferry chlorophyll record was used when tuning the biological model to predict the spring bloom because it is such high frequency, high quality data. It was not possible, however, to tune the model so that each modelled bloom coincided exactly with the ferry record of chlorophyll, and there are slight arrival time discrepancies between the two data sets (Figure 3.20). The 2003 modelled bloom arrives coincident with ferry records in late March, but the ferry shows a sharp decrease in chlorophyll followed by a larger magnitude bloom in late April that the model does not show. The observed second bloom, however, is possibly a different species than the one being modelled, or related to complicated zooplankton dynamics that are not included in this model. The ferry and model agree well in 2004, with the increase in phytoplankton occurring simultaneously in mid-March, although the modelled bloom magnitude is ~10 mg.m~ 3 less than what the ferry recorded. In 2005 the model predicts that the spring bloom arrives about eight days before the ferry records it's arrival, and again, the magnitudes differ. In all three years there is considerable phytoplankton activity following the bloom, including late summer and fall "mini-blooms" in 2003 and 2004, and the model does a moderate job predicting these blooms. 93 Chapter 4 Discussion 4.1 Interannual Variability in the Strait of Georgia The four years sampled during the STRATOGEM period (2002, 2003, 2004, 2005) varied greatly meteorologically. The arrival time of the spring chlorophyll bloom also varied annually, peaking earliest in 2005 (February 26) and latest in 2002 (April 10). The significance of the observed year-to-year differences in wind magnitude, Fraser River flow, and cloud cover, particularly during the spring, are a main focus of this project and necessitate a brief discussion. Data for 2001 is sometimes included because the model uses data from 2001, although that year was not part of the field program. 4.1.1 Wind Mean wind for 2001 through to the end of the field program in June 2005 was ~5 m s - 1 . Wind was strongest in January through March in 2002 (5.62 ms - 1 ) , followed closely by 2003 (5.04 ms - 1 ) and 2004 (5.07 ms - 1 ) , and was significantly weaker in 2005 (3.24 m s - 1 ) (Figure 4.1); the 30-year mean for this period is 3.38 m s - 1 (YVR, 2005). This trend is the same when the period from January 1 to April 1 or May 1 is considered. Wind is clearly strongest in 2002, and weakest in 2005 during the spring; the bloom is latest in 2002 and earliest in 2005 (Figure 4.2). Chapter 4. Discussion 94 10000 - 9000 - 8000 - <f" 7000 - Jan Jul 2003 Jul 04 Jul 05 Jul Figure 4.1: Wind magnitude cubed for the whole sampling period measured at Sand Heads (Figure 1.2). Wind variability is evident, and strongest winds occur in 2002. 4.1.2 Fraser The Fraser River freshet varied greatly over the sampling period (Figure 4.3); the largest freshet occurred in 2002, while the other years' peaks are very similar to one another in magnitude. The 2002 freshet was greater than the climatological mean by more than S x l O W s " 1 (BC Rivers, 2005). Flow in the winter of 2005 is higher than any other year, staying consistently above the winter average of 8 .5xl0 2 m 3 s - 1 , and four times greater than mean winter flow (BC Rivers, 2005) (Figure 4.4 bottom panel). Heavy precipitation in January 2005 is responsible this high flow: January 2005 total precipitation was 249.6 mm compared to a 30-year mean of 153.6 mm. In each year the freshet does not begin until after the bloom is well underway, or over (Figure 4.4). The high flow in winter 2005 cannot be considered a part of the freshet because the flow returned to just above winter levels by the middle of February. Chapter 4. Discussion 95 (a) 2002 (d) 2005 i i . SS! " A 4000 2000 Apr May Figure 4.2: Daily maximum wind magnitude cubed (thin lines) and daily maximum modelled surface chlorophyll (thick lines) for February 1 through. March 30 of each year. Winds are strongest in 2002 and weakest in 2005; the bloom is latest in 2002 and earliest in 2005. Chapter 4. Discussion 96 Figure 4.4: Fraser River flow (thin lines) and modelled daily maximum surface chlorophyll (thick lines) from January 1 to May 1 for each year. Each year the bloom begins before the freshet; the anomalously high flow in January 2005 is apparent. Chapter 4. Discussion 97 3501 1 1 1 1 - — i 1 1 1 1 1 r Figure 4.5: Modelled daily maximum irradiance (IPAR) for 2002-2005, based on cloud cover that is measured at Y V R . Light is high in summer and low in winter; no interannual variations are visible. 4 . 1 . 3 C louds Cloud cover in the Strait of Georgia region typically follows a pattern of obscured skies (cloud cover = 9) throughout much of the winter months and clear skies (cloud cover = 0) during the summer. Average summer (winter) IPAR is ~325 W m ~ 2 (~75 W m - 2 ) (Figure 4.5). February of 2005, the year of the earliest bloom, was much less cloudy than the other years (Figure 4.6), and had 139.9 hours of bright sunshine, a 30-year maximum and far greater than the 30-year mean of 84.6 hours (YVR, 2005), February total precipitation was 45.8 mm compared to a 30-year mean of 123.1 mm (YVR, 2005). Other than this anomaly the four sampled springs were not drastically different from each other in terms of light and cloud cover. Chapter 4. Discussion 98 2 4 0 , | I 1 1 1 1 1 1 r 01/06 01/13 01/20 01/27 02/03 02/10 02/17 02/24 03/03 03/10 03/17 03/24 Figure 4.6: Modelled daily maximum irradiance for January, February and March 2002 and 2005. IpAR is significantly higher in spring 2005 compared to spring 2002. 4.2 What is controlling the spring bloom? To determine the effects of different physical forcings on the arrival time of the spring bloom, two series of tests were run using various data sets. For these tests one physical force was varied and all other variables remained constant. The four model runs (R l , R2, R3 and R4, Table B.2) were run with one variable- either clouds (c), Fraser River (/), or wind ( i t ; ) - changed. The variable is replaced with artificial data or data from another sampled year (Table B.3). The first part of the ID refers to the model run, the second part refers to the variable that is being changed, and the third part refers to the data it has been replaced with. For example, run "RlcR4" means that R l was run with clouds from R4, and everything else was kept the same as in R l ; "R4w0" means R4 was run with winds set to zero; "R2fhalf" means R2 was run with half of R2 Fraser River flow. Tests were also run with Fraser River flow shifted (fshift) ahead by fifty days, such that the freshet peaked much nearer to average bloom peak date. The data sets that were run with data Chapter 4. Discussion 99 from other sampled years were less extreme and more realistic tests of model response to different forcings than the artificial data sets, such as "...c9" meaning clouds were replaced with a constant cloud cover of 9. Results from these tests are shown and discussed in the following sections. Figures 4.8, 4.11, and 4.13 were created by integrating the physical force (wind magnitude cubed, cloud fraction, or Fraser River flow) between January 1 and March 1 for each data series tested. The results of all the tests that were run are plotted against the associated peak bloom date and a regression was done to assess the statistical significance of the relationship. The same figures were created and regressions completed using different time periods, such as January 1 to May 1, February 1 to May 1, February 1 to Apri l 1, and so on. For all three forcings being considered, the time period used produced qualitatively the same results as the figures presented here. 4.2.1 Winds Variations in the wind have a strong effect on the arrival time of the spring bloom, with lower wind producing an earlier bloom and higher wind causing the bloom to be delayed (Figure 4.7). Low wind mixing in the Strait of Georgia provides the stratification necessary for a bloom to develop, whereas strong wind mixing prevents stratification. Between realistic average wind speeds (>1.5.ms - 1 and <6 m s - 1 ) there is a strong linear relationship between wind strength and bloom time arrival (Figure 4.8). The relationship follows the equation Dateof bloom peak — 15.79 days/ms~l x hourly average wind magnitude + Jan.8 (4.1) and it has an R 2 value of 0.75. This result indicates that wind explains 75% of the variance in the spring bloom arrival time. The earliest bloom arrivals, caused by running the model with no wind ("wO"), and the results of the "wl5" model tests are not included in this fit because they Chapter 4. Discussion 100 are unrealistic wind speeds and therefore bias the fit. These results are still valid in a qualitative sense as they are consistent with the outcome of the model testing, i.e.: that wind controls the bloom. Model runs forced with no wind ("wO") caused the bloom to arrive earlier than.any other wind scenario for each year tested, up to more than two months earlier than actual bloom arrival times (Table B.3). The "whalf" tests caused blooms to arrive second earliest each year. Model runs forced with 2005 wind ("wl4, w24, w34") caused the bloom to arrive third earliest- in each case earlier than true arrival time- up to three weeks earlier than actual bloom arrival times. Winds from 2002 ("w21, w31, w41") caused the blooms to occur later than the observed blooms. Winds multiplied by one and a half times the actual wind ("wl5") caused blooms that arrived last three out of the four tested years. Multiplying the wind by two ("...wtwice") resulted in unstable and non-physical/non-sensical results, indicating that there is a point when there is too much wind for the model to accurately process. These findings further enforce the statement that the spring phytoplankton bloom depends on the wind, and will arrive earliest with weakest wind and latest with strongest wind. Once the bloom is in progress wind events can slow the bloom down by mixing phytoplankton down below the critical depth where growth is inhibited by the absence of light (Figure 4.2, March 20, 2004, Figure 3.19). 4.2.2 C louds Variations in cloud cover have a weak effect on the the timing of the spring bloom (Figures 4.9, 4.11). When the model was tested with each different scenario (Table B.3) the bloom peak dates are all within a few days of each other (Figures 4.9). There is, however, a trend evident within each year, suggesting that although wind controls the timing of the bloom, the light does have an Chapter 4. Discussion 101 Figure 4.7: Daily maximum surface chlorophyll (thick line) with daily maximum wind magnitude cubed (thin line) for various wind scenarios between Jan.l and July 1. Arrows indicate the peak of the modelled 2003 bloom with true (2003) wind. No wind (a) and half wind (b) result in the earliest blooms (peaks on Jan.19 and Feb.l, respectively), strongest winds (e) result in the latest bloom (peak on April 1). The bloom in 2005 (d) is much earlier than the 2002 (c) bloom. Chapter 4. Discussion 102 Average hourly wind magnitude (ms* Figure 4.8: Effects of the wind magnitude on the timing of the spring bloom. Date of bloom peak for each wind forcing test is plotted against the hourly average wind magnitude of each forcing test. Bold symbols represent the real data for each year and the corresponding observed bloom time. Chapter 4. Discussion 103 effect. This effect can be quantified by performing a regression on the cloud data with the effect of the wind removed; the result is a prediction of bloom arrival time based on wind and clouds Date of bloom peak — 15.79 days/ms~l x hourly avg. wind magnitude+1.89 x hourly avg.cloud fraction-\-D • (4-2) The fit does a good job predicting the arrival date of the bloom peak compared to the model predicted date of the bloom peak based on known cloud cover and known wind speed (Figure 4.10). Running the model with data from other years caused slight differences in bloom arrival times, as all bloom peaks fell within several days of each other. Light testing results were consistent in producing blooms that followed a logical order, with higher light levels ("cO") producing earlier blooms and lower light levels ("chalf","c9") producing later blooms. "cO" tests caused the bloom to arrive earliest - up to 12 days earlier than actual bloom- in four of the four tested years. Forcing tests "chaff" caused blooms that arrived within 2 days of the observed bloom each year, and total cloud cover ("c9") caused blooms to arrive between 1 and three days later than observed blooms in all but one year. These results are evidence of the correlation between cloud cover and peak bloom date described above. 4.2.3 Fraser River Variations in Fraser River flow do not effect a consistent pattern in bloom arrival time (Figures 4.12, 4.13). Tests were run using multiples of Fraser River flow, as well as a shifted flow so that the freshet arrived fifty days earlier- much closer to a typical bloom time arrival. There is a strong scattered nature to Figure 4.13, showing that no statistically significant relationship between flow and bloom arrival exists. None of the tests in Table B.3 produced any conclusive results (Figure Chapter 4. Discussion 104 500 Figure 4.9: Daily maximum surface chlorophyll (thick line) with modelled daily maximum irradi- ance (thin line) for various cloud cover scenarios between Jan.l and July 1. Arrows indicate the peak of the 2004 bloom with true (2004) cloud cover. The bloom arrives a few days earlier with clear skies (a) than for completely obscured skies (d), showing that light intensity has a very small (insignificant) effect on the bloom arrival time. The sunny skies in February 2005 (e) do not cause an early bloom, implying that the bloom arrival time does not depend on light except in extreme cases (constant cloud cover, e.g. (a)). Chapter 4. Discussion 105 Date predicted by fit for clouds and wind Figure 4.10: Comparison of bloom arrival date predicted by (4.2) to dates predicted from model testing with different cloud scenarios. 4.12), therefore Fraser River flow has no effect on bloom arrival time. 4.3 Tides The model includes tidal mixing implicitly in the internal parameterization as part of the shear diffusion. The tides in the Strait of Georgia are important to the surface salinity due to the movement of the plume by the tide (Stronach, 1977). The ferry data reveals that the tidally advected plume crosses over S3 with the ebb and flood tides and changes the surface salinity there (Figure 3.13). As tidal effects are included only implicitly, it is important to verify that the tidal influence on salinity is modelled effectively. To see if the strong effects of the tide on the plume were being missed by the model the tidal excursion (the path of an average water particle during Chapter 4. Discussion 106 4 5 6 Average hourly cloud fraction Figure 4.11: Effects of cloud cover on the timing of the spring bloom. Date of bloom peak for each cloud forcing test is plotted against the hourly average cloud fraction of each forcing test. Bold symbols represent the real data for each year and the corresponding observed bloom time. Chapter 4. Discussion 107 Figure 4.12: Daily maximum surface chlorophyll (thick line) with Fraser River flow (thin line) for various flow scenarios between Jan.l and July 1. Arrows indicate the peak of the 2003 bloom with true (2003) Fraser River flow. There is no variation in the 2003 bloom time arrival with shifted flow (b) or twice the actual flow (e). The bloom is a few days late with half the actual flow (a) and the 2002 flow (d) and a few days early with 2005 flow (c). No pattern is discernible, hence, Fraser River flow has no effect on bloom arrival time. Chapter 4. Discussion 108 May Apr CO °> CL + E o o Mar o B to Q Feb Jan AO © O • * * A A * o 0 2002 A 2003 * 2004 + 2005 500 1000 1500 2000 2500 3000 3500 Average daily flow (m 3 s 1 ) 4000 4500 Figure 4.13: Effects of Fraser River flow on the timing of the spring bloom. Date of bloom peak for each Fraser River forcing test is plotted against the daily average flow of each forcing test. Bold symbols represent the real data for each year and the corresponding observed bloom time. Chapter 4. Discussion 109 the tidal cycle) for S3 was calculated using T.E. — {tirti€.iQWWateT ti7ne{ligilwater) x v^ax 2/TT (4.3) where vmax is the maximum tidal current velocity. The value used for vmax is the maximum tidal current velocity at SI (Figure 2.1) for the M2 tidal component semi-major axis at the surface (Foreman et al., 1995). The mixed semi-diurnal nature of the tide in the Strait of Georgia (Figure 1.6) means that the tide goes from HHW to HLW to LHW to LLW and back to HHW, so two values for t«me; o t u l U Q ter ~ timehighwater f ° r the ebb tidal excursion and two values for t i m e i o w w a t e r — timehighwater on the flood tidal excursion were calculated. The result is that a particle can travel up to 7 km towards the Pacific on large ebb, and 7 km up-Strait on a large flood, and ~4 km on a small ebb or flood. This translates into an along-latitude variation (assuming no longitudinal movement occurs) of 0°03.43 (~4 minutes of latitude). To compare the salinity of a particle travelling on the tide to particles that are stationary at S3, as in the model, we can plot ferry salinity data from S3 (one minute on either side of S3) and ferry salinity data within the tidal excursion (4 minutes on either side of the latitude, but within the longitude of S3) (Figure 4.14). We see that the salinity can vary by as much as 10 on a flood tide and 5 on an ebb tide. The implications of this dependence of surface salinity on the tide are that the salinity recorded when a CTD cast is taken on a monthly survey will be a reflection of the tidal cycle during which the data was recorded. The model does not account for these variations "in surface salinity, so when comparing modelled surface salinity to observed surface salinity in profiles like Figure 3.3 we need to consider that if the CTD cast was taken on a flood tide the surface salinity could be up to 10 greater than what the model computes for surface salinity. Analysis of tidal records and Chapter 4. Discussion 110 Figure 4.14: The ferry recorded salinity within the calculated tidal excursion box and ferry recorded salinity at S3 showing the variability of the salinity due to tidal fluctua- tions. surface salinity discrepancies between modelled and observed surface salinity revealed that, indeed, some of the occasions when observations were more saline than the model can be explained by a flood tide, and some of the occasions when observations were much fresher than the model can be explained by an ebb tide. However, the converse was also found to be true such that salinity was under-predicted during a flood tide and over-predicting during an ebb. These findings are not necessarily contradictory, as we know the Fraser River plume is very spatially variable both with the tide and the wind (Halverson et a l , 2005; Hoos and Packman, 1974; Stronach, 1977). Ultimately these variations in salinity are on very small timescales and will not affect the timing of the spring phytoplankton bloom- the focus of this thesis. Chapter 4. Discussion 111 4 . 4 Model Performance Summary 4.4.1 M o d e l L imi t a t i ons As with any modelling study many of the results of this research are limited.by the constraints imposed by the model itself. Although the model was optimized for the Strait of Georgia as much as possible it was also simplified as much as possible, and these simplifications put significant limitations on the accuracy and results. The two major limitations to consider are the geographical separation of the data sources and their distance from the actual modelled location, and the re- parameterization of the physical model. The physical model was limited by spatial constraints. Using a 2 or 3-Dimensional model introduces much more complexity but also represents a scenario that is more similar to the real physical environment being modelled. Spatial variability in the Strait is lost by modelling only one finite point (Station S3) rather than a series of stations, but in this case the benefits of using a simpler model outweigh the loss of complexity. The clear results of this study might not have emerged with a more complex model, so the spatial limitation is an acceptable one. Each re-parameterization done to optimize the physical model for the Strait of Georgia has within it some limitations. The freshwater flux and light parameterizations were both extensive and often based on innovative methods and a wide range of data. The surface salinity in the model is calculated based on a regression of daily flow estimates and matched to the fifteen minute model time step in a smoothing procedure. This value is used to calculate the freshwater flux using an equation that is based loosely on the Knudson relations. Al l these steps introduce the possibility of error and inaccuracy. Another issue related to the freshwater flux is the Fraser River flow data, which is of particular Chapter 4. Discussion 112 concern as it is measured more than 1x10 s m away from S3. However, the spring bloom is not sensitive to the Fraser River flow, so this is not overly problematic. The light parameterizations are of greater importance to the bloom dynamics because of the phytoplankton's dependence on light and heat. The cloud model was completely altered using cloud data measured over land ( Y V R ) and comparing it to irradiance data measured over land several kilometres away (UBC) in an attempt to predict light over water despite known differences between clouds over land and clouds over water. For example, the land-water gradient on a warm summer afternoon often causes clouds to form over land, whereas the sky over the ocean is cloud free. The albedo used in the model was calculated using crude equipment and an experimental procedure, but we can see from the comparison with the Halibut Buoy data that the parameterization of surface IPAR including the clouds has been successful. The parameterization for light below the water is somewhat more uncertain. The calculation for P A R is based on four variables, two of which are calculated by the model and therefore have all the error and uncertainty inherent in those calculations. The model chooses how total light decreases with depth according to Jerlov (1976) based on a best-fit algorithm, although the Jerlov data itself is a best-fit generalization. Irradiance profiles are often under-predicted (light gets too deep), an error likely due to the generally under-predicted P that K depends on (2.30). Model results for temperature are good even given the limitations to the total light parameterization, however, so it would appear that ultimate model performance is not constrained by any oversights and generalizations inherent in that adaptation. The parameterization of the pressure gradients is approximate, and depends on several assump- tions including the central location of station S3, nearly vertical sides of the basin, homogeneity, etc. However, what matters in the model is the strength and depth of the currents, not the direc- Chapter 4. Discussion 113 tion, and the parameterization likely does a much better job of the former. Therefore the pressure gradient scheme does not limit the model performance significantly. The biological model was simplified greatly by including only one source of nutrients and no regeneration, one group of phytoplankton, and no zooplankton. It is known that other groups and sizes of phytoplankton are important in the ecology of the Strait of Georgia (Parsons et al., 1978). Flagellates are especially important, dominating the phytoplankton community during the winter. The lower-than-observed magnitude of modelled phytoplankton could be explained in part due to the absence of the interaction between the different size classes in the Strait of Georgia. Zooplankton, especially copepods, also play an integral role in bloom dynamics in the Strait of Georgia (Harrison et al., 1983) that, aside from a background loss term that includes respiration and mortality, are largely ignored in the model. The timing of the arrival of the copepods at the surface, their abundance, their grazing rate- these are all factors which play an important role in the arrival, magnitude, and termination of the phytoplankton bloom. Typically, the copepods do not arrive until after the arrival of the spring bloom, so the model does a good job predicting this arrival despite the lack of a dynamic zooplankton component. However, improvements to modelled bloom magnitude and phytoplankton levels during the rest of the year should be much improved by the addition of an active zooplankton component. The absence of ammonium and therefore the regeneration loop is also a source of uncertainty. As nitrate is plentiful leading up to the bloom, this process is inconsequential, and therefore the model does a good job predicting the timing of the spring bloom. The contribution of sinking biomass to the regeneration loop would add to the phytoplankton food sources and result in greater growth, leading to larger modelled blooms. Another major limitation of the nutrient modelling is the utilization of a half-saturation constant for nitrate in the nutrient uptake equation. Instead, a Chapter 4. Discussion 114 half-saturation uptake for growth should have been used. This oversight can perhaps account for some of the low modelled nitrate levels. The implications of leaving these complexities out of the biological model are that the ecosystem is not represented as accurately as possible, and crucial information could be lacking in the modelled system. The disagreement between the modelled magnitude of the spring bloom and the observed magnitude can be accounted for by these limitations. Inclusion of different types of nutrients, other species, and a more complex grazing scheme could bring the chlorophyll concentrations up to a level that agrees with the observed level and concurrently improve the IPAR parameterization. The non-temperature dependent nature of the loss term, RM, is also a source for uncertainty, as the inclusion of such a process would add seasonality to this term. As the model is most sensitive to this parameter, a temperature dependency would be an important process to include in order to correctly predict the bloom timing and magnitude. 4.4.2 Cons is tency of Resul ts The physical model predicts vertical profiles of temperature, salinity, density, and irradiance at Station S3 in the Strait of Georgia that are consistent with observed profiles of these variables. Mixed layer depth, pycnocline, light attenuation, eddy viscosity and diffusivity are well-modelled. The physical model responds appropriately to external forcings, such that an increase in wind speed deepens the mixed layer depth (Figure 3.7) and increases surface currents, and increased freshwater input causes the mixed layer depth to become shallower and the surface salinity to become fresher (Figures 3.7, 3.8). The coupled physical and biological model predicts the spring phytoplankton bloom within a few days of the observed bloom (Figure 3.16) and the termination of the bloom coincides with the Chapter 4. Discussion 115 modelled and observed nitrate depletion (Figure 3.17). The biological parameters also respond appropriately to external forcing such that increased wind and the associated deepened mixed layer causes the plankton to be mixed out of the surface layer to deeper depths (Figure 3.18); when P are mixed below the critical depth there is a decrease in biomass (Figure 3.19). 4.4.3 C o m p a r i s o n w i t h Other Resul ts Another way to evaluate the success of the model and significance to the scientific community is to compare it to similar contemporary models or studies that strive to achieve the same goals. A comparison of the results of this study to similar models or studies is also a good indication of the success of the study. Based on observed biological and physical data Yin et al. (1996) asserted that the bloom arrival is dependent on the stratification produced by the Fraser River freshet. This theory has been strongly refuted by the results presented in this thesis, which clearly show that the Fraser River flow has no effect on the bloom' arrival time (Figure 4.13, 4.12). A bloom will occur when the water column is stratified, and I assert that this stratification is achieved by the absence of wind mixing, contrary to Yin et al. (1997b), who claim that the necessary stratification for a bloom can only be produced by river discharge. If wind is weak, state Yin et al. (1997b), even a low river discharge can produce stratification, but if wind is strong, the bloom will not occur until river discharge overcomes wind mixing. Using a high-resolution mixing model developed for the Strait of Georgia to perform extensive model testing, I have shown that, although the empirical evidence may suggest an important role for the Fraser River in creating stratification, the bloom will occur only when winds are weak, regardless of the river discharge. The relationship between wind and bloom time arrival discussed in this thesis agrees .with Yin's claim that wind can delay the bloom Chapter 4. Discussion 116 (Yin et al., 1997b). The coupled physical and biological two-layer box model of L i et al. (2000) found that plankton populations were insensitive to changes in estuarine circulation caused by variability in Fraser River discharge, a result that correlates well with the results of this thesis. L i et al did not include wind in their model, nor did they test the model's sensitivity to changes in light availability. The conclusion of Jin et al. (2005) that the absence of wind caused the bloom in the Bering Sea to be delayed by several weeks does not agree with the conclusion of this study, that no wind causes an earlier bloom. This disagreement of results could stem from the fact that their model was based on only one year of data, whereas the Strait of Georgia model predicted the bloom for four years, or on the different environmental setting of their model. Predicting the bloom under four meteorologically and biologically different years adds certainty to the model that wouldn't exist for a one year prediction. The model of Eslinger et al. (2001), which predicts five spring blooms within weeks of observa- tions, ran from late February to mid-November for 1986 through 1998. This model differed from the Strait of Georgia model primarily in that all simulations began with identical N,P,Z and T homogeneously distributed throughout the water column. This was done to isolate the ecosystem from the physical environment to determine the extent to which meteorological forcing explained spring bloom variability. Their model also included neither salinity nor advection, both key pa- rameterizations in the Strait of Georgia model. Solar radiation (IPAR and ITOTAL) w a s calculated by Eslinger et al. (2001) using a model with cloud concentration set constant at 35%. Although the study compared several years of extensive observations the authors did not use their model as a tool to test their conclusions as we did. Their results are based solely on the strength of the constant initial profile parameterization and on observations. Chapter 4. Discussion 117 The linear and stepwise regression analyses of Iriarte and Purdie (2004) found no clear rela- tionship between wind speed and the timing of the spring bloom. The study considered daily averaged wind speed for three days prior to the bloom, which was defined as the first time the chl a concentration was higher than 10 mg m~ 3 . This study is a pure observational analysis and does not compare to the model used here, although the inclusion of several different years and a complete set of data adds value to the conclusions. 4.5 Implications 4.5.1 Climate There are many factors that play a role in the complicated climatology of the wind system in the North Pacific Ocean, from large-scale ocean-atmosphere systems to climate change. If global temperatures rise as predicted scientists expect a reduction in global wind fields due to decreased temperature gradients from low to higher latitudes (National Coastal Assessment Group, 2000). However, the possibility exists that land-ocean temperature gradients will change with increased heating over land and cause increases in local wind (National Coastal Assessment Group, 2000). These, scenarios are conflicting with respect to their implications for the ecology of the Strait of Georgia, with the former potentially causing earlier blooms and the latter potentially causing later blooms. The Aleutian Low Pressure System controls atmospheric circulation in the North-East Pacific between September and May and is responsible for many of the storm systems in that region (Overland et al., 1999). The intensity of the mechanical forcing of the ocean can be determined from the strength of the Aleutian low by analysis of the surface pressure gradient (Overland et al., Chapter 4. Discussion 118 1999). The system has undergone several regime shifts on a decadal scale over the past century; these are explained by variability in other atmosphere-ocean systems (primarily the Pacific-North America (PNA) pattern and the Arctic Oscillation). The variations in the strength of the PNA cause variations in the Aleutian low strength and position, which is linked to wind, temperature and precipitation patterns (Stabeno et al., 2001). Predictions for trends in the P N A are conflicting, with some models predicting an increase and others predicting a decrease in the P N A over the next century (National Coastal Assessment Group, 2000). Most regional climate model grids are too coarse to resolve the small-scale winds in the Strait of Georgia, but some general trends for the Northeastern Pacific can be extracted. The first simulation of the Canadian Regional Climate Model (Laprise et al., 1998) predicts a drop in winter Mean Sea Level Pressure over the Strait of Georgia of 5 hPa increasing to 7 hPa off the west coast of Vancouver Island when C O 2 concentrations are doubled. This indicates a deepening and southward displacement of the Aleutian Low, which suggests stronger winds for the Strait of Georgia. The most recent version of the C R C M Canada (2005) suggests a weakening of meridional winds and a strengthening of zonal winds between 1994 and 2063 in the Northeastern Pacific on the order of 1 m s - 1 under an enhanced greenhouse gas scenario. Changes over the Strait of Georgia are insignificant. Wind trends over the STRATOGEM sampling period suggest that wind strength is increasing, as 2002, 2003 and 2004 are all higher than the 30-year mean wind speed by more than 1.5ms -1. The effects of climatic shifts such as E l Nino and La Nina include changes in wintertime air temperature, precipitation, alongshore wind, and sea level; these tend to increase for the former and decrease for the latter system (Stabeno et al., 2001). The E l Nino system is currently classified as neutral but there is considerable uncertainty surrounding this system (Climate Prediction Center, Chapter 4. Discussion 119 2005). Thus, there are many factors that combine to create the wind system in the Strait of Georgia and a prediction of what may happen would be a gamble at best. The effect that an earlier spring bloom has on the rest of the ecosystem is the next big problem facing the S T R A T O G E M group. 4.6 Future Work The full significance of this research cannot be realized without further work, both on the model and in other areas of the project. These include putting back some of the variables that were taken out of the biological model, such as other nutrient sources, other phytoplankton species, and a zooplankton component. The addition of any one of these will increase the validity of the results of this thesis because it will then be a closer representation of the system being modelled. It is also hoped that these additions will improve model performance, specifically in terms of the magnitude of the bloom and the fall nitrate replenishment. Much of the other work that remains to be completed on this topic is currently underway as other components of the STRATOGEM project. Of particular interest to this research is the implications of the timing of the spring bloom for the upper trophic levels, from zooplankton up to commercial fisheries. The dynamics of zooplankton in the Strait of Georgia are being investigated by the University of Victoria component of the S T R A T O G E M project; those results will enable further, more comprehensive studies on the implications of lower trophic level dynamics and physical environmental effects on the whole ecology of the Strait of Georgia. 120 Chapter 5 Conclusion The objective of this thesis was to elucidate the relationship between physical forcings and the timing of the spring phytoplankton bloom, particularly the effects of wind, sunlight and freshwater flux (primarily from the Fraser River) on the arrival date of the bloom. This was a question that had been posed to the scientific community before but conclusive evidence was sparse. A "set of rules" was sought to govern the biological and physical interactions of the Strait of Georgia. The conclusions of this thesis are meant to be the first set of rules, dictating what determines the arrival time of the spring bloom, so that ultimately one can make assumptions on ecological responses of the Strait based on climatological indicators. A ID vertical mixing model was adapted from a well-known K P P model (Large et al., 1994) and optimized to model the Strait of Georgia. Surface salinity and freshwater flux parameterizations are based on Fraser and Englishman River fluxes; vertical advection is also controlled by the Fraser River. The enclosed nature of the Strait of Georgia requires the addition to the model of a baroclinic pressure gradient parameterization as the pressure-driven currents contribute to the overall circulation of the estuary. Light above the water is parameterized using a localized cloud model and an experimentally determined albedo value. Light below the water is separated into the total light spectrum which contributes to the non-turbulent heat flux parameterization, and Photosynthetically Active Radiation which is required for phytoplankton growth. Total light is Chapter 5. Conclusion 121 absorbed in the water column following Jerlov (1976) and a dependence on Fraser River flow and phytoplankton concentration, while biologically available light is absorbed based on an empirical fit of the attenuation coefficient (K) to Fraser River flow, phytoplankton concentration, and the mixed layer depth. Modelled results of salinity and light compare well to observational data indicating that the processes governing each has been parameterized successfully. The physical model was then coupled to a simple biological model that included nitrate as the only nutrient and microphytoplankton as the only plankton. The phytoplankton being modelled is the most dominant phytoplankton in the Strait, the chain-forming diatom Skeletonema costatum, which is the main species in the ~30 mg m - 3 bloom of each spring. Ammonium and urea are less readily available in the Strait than nitrate, and so their exclusion from the biological model is justified. Zooplankton are included only as a grazing constant. The model was tuned to correctly predict the arrival of the spring phytoplankton bloom by adjusting values of losses to phytoplankton growth (Rm). The model is initialized with profiles of temperature, salinity, nitrate, and fluorescence mea- sured as part of the S T R A T O G E M monthly field sampling program. Hourly wind, air temperature, cloud fraction, and humidity data force the model in addition to daily Fraser and Englishman River flow. Temperature, salinity and fluorescence data measured on a B C Ferries vessel as part of the S T R A T O G E M ferry sampling program is used for model comparisons. The model correctly predicts the arrival time of the bloom in each of 2002-2005 and the corresponding associated nitrate depletion. The modelled bloom is of lesser magnitude than ob- servations due to the absence of a more complex biological model that includes additional nutrient sources, phytoplankton species, and an active zooplankton component. Modelled nitrate matches observations well throughout the winter and spring but surface nitrate does not increase in the Chapter 5. Conclusion 122 fall as it should due in part to the existence of stronger than observed stratification. Profiles of physical and biological parameters compare reasonably well to observations. The mixed layer depth is controlled throughout the winter by wind and is frequently deeper than 5 m, but the stratification produced by the Fraser River freshet maintains a shallow mixed layer depth (less than 5 m) throughout the late spring and summer. The model was forced with varying wind, cloud cover and freshwater flux data sets to identify and investigate any relationship between physical environmental force and biological response. It was found that wind controls the arrival time of the spring phytoplankton bloom in the Strait of Georgia, explaining 75% of the variance in bloom arrival time; cloud cover has a secondary effect on the timing of the bloom. The date of the peak of the bloom can be predicted based on average hourly wind speed and hourly average cloud fraction for January and February of a given year. Freshwater input from the Fraser River has no effect on the timing of the bloom whatsoever. Future climatic shifts will undoubtedly have some effect on the wind system of the Strait of Georgia, which will then effect a change in the spring bloom dynamics there. The implications of an earlier or later than average bloom on the upper trophic levels remains to be discovered, and will surely continue to challenge the oceanographic community long after S T R A T O G E M is complete. 123 Bibliography Ascher, U . , S. Ruth, and B. Wetton, 1995: Implicit-explicit methods for time-dependent partial differential equations. SIAM Journal of Numerical Analysis, 32, 797-823. Batchelder, H . , 2003: U.S. Globec: Who we are and what we do. World Wide Web, http : I/www .usglobec.org/gaag/overview .html. B C Rivers, 2005: Real-time hydrometric data. World Wide Web, http : //scitech.pyr.ee. gc.ca/waterweb/main.aspllang = 0. Beamish, R. J . , C . E . M . Neville, and A . J . Cass, 1997: Production of Fraser River sockeye salmon (Oncorhynchus nerka) in relation to decadal-scale changes in the climate and the ocean. Can. J. Fish.' Aquat. Sci., 54, 543-554. Beaugrand, G . , P. C . Reid, F . Ibanez, J . A . Lindley, and M . Edwards, 2002: Reorganization of North Atlantic marine copepod biodiversity and climate. Science, 296, 1692-1694. Biospherical Instruments Inc., 2003: Frequently Answered Questions. World Wide Web, http : I/www. biospherical.com/. Boom, B. D. L . , 1976: Mathematical modelling of the chlorophyll distribution in the Fraser River plume, BC. Master's thesis, U B C . Bibliography 124 Canada, E . , 2005: Canadian Center for Climate Modelling and Analysis. C R C M 3 . 6 Data. World Wide Web, http : //'www.cccma.bc.ec.gc.ca/'data/'crcm36/crcm36.shtml. Canadian Heritage Rivers System, 2005: Fraser River, British Columbia: Where the Salmon is King. World Wide Web, http : //www.chrs.ca/Rivers/Fraser/Fraser — Fe:htm. Cassis, D. , 2005: Notes on S T R A T O G E M spring blooms 2004 and 2005, data observations and planktonic analysis. Chen, C , 2004: Marine ecosystems dynamics modeling. World Wide Web, http : //codfish.smast.umassd.edu/. City of Parksville, 2005: City of Parksville. Operations - Public Works. World Wide Web, http : //city, parksville.be. ca/. Climate Prediction Center, 2005: E l Nino Southern Oscilla- tion (ENSO) Diagnostic Discussion. World Wide Web, http : //www .cpc.ncep.noaa.gov/products/analy sismonitoring/ensoadvisory/. Cloern, J . E . , 1978: Empirical model of Skeletonema costatum photosynthetic rate, with appli- cations in the San Francisco Bay Estuary. Advances in Water Resources, 1, 267-274. Dobson, F . W . and S. D. Smith, 1988: Bulk models of solar radiation at sea. Q. J. R. Meteoro- logical Society, 114, 165-182. Edwards, M . and A . J . Richardson, 2004: Impact of climate change on marine pelagic phenology and trophic mismatch. Nature, 430, 881-883. Bibliography 125 Eppley, R. W . , 1972: Temperature and phytoplankton growth in the sea. Fishery Bulletin, 70, 1063-1081. Eslinger, D. , R. Cooney, C . McRoy, A . Ward, T . K . Jr., E . Simpson, J.Wang, and J . Allen, 2001: Plankton dynamics: observed and modelled response to physical conditions in Prince William Sound, Alaska. Fisheries Oceanography, 10, 81-96. Falkowski, P. G . , 1975: Nitrate uptake in marine phytoplankton: Comparison of half-saturation constants from seven species. Limnology and Oceanography, 20, 412-417. Fisheries and Oceans, Canada, 2005: Fisheries Management - Pacific Region. World Wide Web, http : //www.pac.df o — mpo.gc.ca/ops/fm/fishmgmte.htm. Foreman, M . , R. Walters, R. Henry, C . Kellers, and A . Dolling, 1995: A tidal model for Eastern Juan de Fuca Strait and the Southern Strait of Georgia. Journal of Geophysical Research, 100, 721-740. Gargett, A . , 1972: Proceedings of the Strait of Georgia Workshop. Technical report, Institute of Oceanography, University of British Columbia. Gower, J . , 2005: Data from Halibut Bank Buoy, courtesy of Fisheries and Oceans, Canada. Halverson, M . , R. Pawlowicz, R. Lee, and S. Allen, 2005: High resolution ferry-based observa- tions of three consecutive spring blooms in the Strait of Georgia. Canadian Meteorological and Oceanographic Society. Harrison, P., J . Fulton, F . Taylor, and T . Parsons, 1983: Review of the biological oceanography of the Strait of Georgia. Can. J. Fish. Aquat. Sci., 40, 1064-1094. Bibliography 126 Hobson, L. A. and M . R. McQuoid, 1997: Temporal variations among planktonic assemblages in a turbulent environment of the Southern Strait of Georgia, British Columbia, Canada. Marine Ecology Progress Series, 150 , 263-274. Hoos, L. M . and G. A. Packman, 1974: The Fraser River Estuary Status of Knowledge to 1974. Special Estuary Series 1, Environment Canada, Vancouver, BC. Iriarte, A. and D. A. Purdie, 2004: Factors controlling the timing of major spring bloom events in an U K south coast estuary. Estuarine, Coastal and Shelf Science, 6 1 , 679-690. Jeffery, N. , 2002: Modelling a Phytoplankton Dichotomy in the Eastern Subarctic Pacific: Impact of Atmospheric Variability, Iron Surface Flux, and Life Cycle Dynamics of the Calanoid Copepods, Neocalanus spp.. Ph.D. thesis, UBC. Jerlov, N . G., 1976: Marine Optics, volume 14 of Elsevier Oceanography Series. Elsevier Scientific Publishing Company, New York. Jin, M . , C. Deal, J. Wang, N . Tanaka, and M . Ikeda, 2005: Vertical mixing effects on the phyto- plankton bloom in the southeastern bering sea mid-shelf, submitted to J. Geophysical Research. Jorgensen, E. G., 1966: Photosynthetic Activity during the Life Cycle of Synchronous Skele- tonema Cells. Physiologia Plantarum, 19 , 789-799. Ketler, R., 2005: Solar radiation data, UBC Department of Agroecology, obtained through per- sonal contact: rick.ketler@ubc.ca. Kobari, T. and T. Ikeda, 2001: Ontogenetic Vertical Migration and Life cycle of Neocalanus Plumchrus (Crustacea: Copepoda) in the Oyashio Region, with Notes on Regional Variations in Body Sizes. Journal of Plankton Research, 23, 287-302. Bibliography 127 Laprise, R., D. Caya, M . Giguere, G. Bergeron, H. Cote, and J. Blanchet, 1998: Climate and Climate Change in Western Canada as Simulated by the Canadian Regional Climate Model. Atmosphere-Ocean, 36, 119-167. Large, W., J . McWilliams, and S. Doney, 1994: Oceanic vertical mixing: A review and a model with a nonlocal boundary layer parameterization. Reviews of Geophysics, 32, 363-403. Large, W. and S. Pond, 1982: Sensible and latent heat flux measurements over the ocean. Journal of Physical Oceanography, 12, 464-482. LeBlond, P., 1983: The Strait of Georgia: functional anatomy of a coastal sea. Can. J. Fish. Aquat. Sci., 40, 1033-1063. LeBlond, P., K . Dyck, K. Perry, and D. Cumming, 1983: Runoff and precipitation time series for the coasts of British Columbia and Washington State. Manuscript Report 39, Department of Oceanography: University of British Columbia. L i , M . , K . Denman, and A. Gargett, 2000: What Determines Seasonal and Interannual Vari- abiulity of Phytoplankton and Zooplankton in Strongly Estuarine Systems? Application to the semi-enclosed estuary of Strait of Georgia and Juan de Fuca Strait. Estuarine, Coastal and Shelf Science, 50, 467-488. L i , M . , A. Gargett, and K. Denman, 1999: Seasonal and Interannual Variability of Estuarine Circulation in a Box Model of the Strait of Georgia and the Juan de Fuca Strait. Atmosphere- Ocean, 37, 1-19. Linguanti, J., 2001: CTD data, obtained through personal contact: LinguantiJ@pac.dfo- mpo.gc.ca. Bibliography 128 Mackas, D. L . and P. J . Harrison, 1997: Nitrogenous Nutrient Sources and Sinks in the Juan de Fuca Strait/Strait of Georgia/Puget Sound Estuarine System: Assessing the Potential for Eutrophication. Estuarine, Coastal and Shelf Science, 44, 1-21. Mount Arrowsmith Biosphere Reserve, 2001: Mount Arrowsmith Biosphere Reserve, South-Eastern Vancouver Island, British Columbia, Canada. World Wide Web, http : //www .mountarrowsmithbiosphere.ca/characteristics .htm. National Coastal Assessment Group, 2000: C O A S T A L : The Potential Consequesnces of Climate Variability and Change. N O A A Coastal Ocean Program, Decision Analysis Series No.21. Overland, J . , J . Adams, and N. Bond, 1999: Decadal variability of the Aleutian Low and its relation to high latitude circulation. J. Climate, 12, 1542-1548. Parsons, T . , 1979: Some ecological, experimental and evolutionary aspects of the upwelling ecosys- tem. South African Journal of Science, 75, 536-540. Parsons, T . , P. Harrison, and R. Waters, 1978: A n experimental simulation of changes in diatom and flagellate blooms. J. exp. mar. Bioo. EcoL, 32, 285-294. Parsons, T . and T . Kessler, 1987: A n ecosystem model for the assessment of plankton production in relation to the survival of young fish. Journal of Plankton Research, 9, 125-137. Parsons, T . R., M . Takahashi, and B. Hargrave, 1984: Biological Oceanographic Processes. Perg- amon Press, New York, 3rd edition. Price, N . , W . Cochlan, and P. Harrison, 1985: Time course of uptake of inorganic and organic nitrogen by phytoplankton in the Strait of Georgia: comparison of frontal and stratified commu- nities. Marine Ecology Progress Series, 27, 39-53. Bibliography 129 Richardson, A . J . and D. S. Shoeman, 2004: Climate impact on plankton ecosystems in the Northeast Atlantic. Science, 305, 1609-1612. Riche, O. , 2005: S T R A T O G E M 2002-2005: Transports and Biochemical Tracer Fluxes Derived with a Three-Box Inverse Model of the Strait of Georgia. Canadian Meteorological and Oceano- graphic Society. Stabeno, P., N . A . Bond, A . J . Hermann, C . W . Mordy, and J . E . Overland, 2001: Meteorology and Oceanography of the Northern Gulf of Alaska, for submission to Continental Shelf Research. Stronach, J . , 1977: Observational and Modelling Studies of the Fraser River Plume. Phd, Uni- versity of British Columbia. Swain, L . , 1998: Water quality assessment and objectives for the Fraser River from Hope to Sturgeon and Roberts Banks. Technical report, Water Management Branch, E n - vironment and Resource Division, Ministry of Environment, Lands and Parks., http : //wlapwww .gov .bc.ca / wat / wq / objectives j hopebank/lower f raser.html. Thomson, R. E . and I. V . Fine, 2003: Estimating mixed layer depth from oceanic profile data. American Meteorological Society, 20, 319-329. Waldichuck, M . , 1957: Physical oceanography of the strait of georgia, british Columbia. J. Fish. Res. Board Can., 14, 321-486. Yin , K . , R. H . Goldblatt, P. J . Harrison, M . A . S. John, P. J . Clifford, and R. J . Beamish, 1997a: Importance of wind and river discharge in influencing nutrient dynamics and phytoplankton production in summer in the central Strait of Georgia. Marine Ecology Progress Series, 161, 173-183. Bibliography 130 Yin, K. , P. J. Harrison, R. H. Goldblatt, and R. J. Beamish, 1996: Spring bloom in the central Strait of Georgia: interactions of river discharge, winds and grazing. Marine Ecology Progress Series, 138 , 255-263. Yin, K. , P. J. Harrison, R. H. Goldblatt, M . A. S. John, , and R. J. Beamish, 1997b: Factors controlling the timing of the spring bloom in the Strait of Georgia estuary, British Columbia, Canada. Y V R , 2005: BC Climate Database. World Wide Web, http : //www.weather of fice.pyr.ec.gc.ca/Climate/. 131 Appendix A Physical Model Wind Stress Wind stress is calculated in the following way Large and Pond (1982) CD = 1.0 x 10- 3( C ° + Cp + C^U20 + V20) (A.1) \Juw + vw where CQ=2.70ms _ 1 , C^O.142 and C 7=0.0764m _ 1s, (Uio.Vio) are the components of the wind velocity at 10 m above the ground. ru = , U W CDPatm\yJu^+V?0\\rv = V w CDPatm\yJu?0 + V?0\2 (A.2) Entrainment Calculations Assuming steady flow, split between layers given by model calculated mixed layer depth, mass balance implies Q 1 = i ? + Q D l (A.3) where R is the river flow and Q i and Qjr>i are the upper and lower layer flows at a given point (Figure A . l ) and therefore the salt balance of the flows gives, since the salinity of the river is 0 QiSi = SDQDl (A.4) Appendix A. Physical Model 132 where Si and So are the salinities of the upper and lower layers. Similarly at a second point Q 2 = R + QD2 (A.5) Q2S2 = SDQD2 (A.6) the upper layer flow at S2 is due to the flow at SI plus any water that has been entrained from below Q2 = Qi + E (A.7) where E is the entrainment flux. These have salt balance S2Q2 = SlQl + SDE (A.8) Eliminating QD1 from (A.3),(A.4) gives Qi = - ^ r (A.9) and similarly for Q2 (On = So - S2 Q2 = ^ r - (A.10) Substituting (A.9),(A.10) into (A.7), we get RSD{S2-SX) ( A n ) Recalling the empirical relation (2.19) but excluding the Englishman River, the surface salinities (Si,S2) can be estimated as S l = SD-- (A. 12) 71 and S2 = SD-- (A. 13) 72 Appendix A. Physical Model 133 for two empirical constants 71 and 72, so that S z > - S ! = - (A.14) 7i S D - S 2 = —. (A.15) 72 Inserting these into (A.11) allows us to eliminate R, and the entrainment is £ = S D ( 7 i - 72). ( A - 1 6 ) Entrainment is upwelling over the area between S2-2 and S3 (Figure A.2), i.e.: g = t , T r ( r V r 2 l ) - (A.17) where u; is the entrainment (upwelling) velocity (2.23) and 1/Cw is the fraction of the area over which upwelling occurs. So a calculation for w leads to the final equation for upwelling velocity, with units of m s - 1 : CwSpj-n - 7 2 ) IR\ w = i~2 21—• (A.18) This is the maximum upwelling velocity and it occurs at 40 m. This value linearly decreases to the surface so that we end up with a value w(l) for each grid level I. Using these values the concentration change of a property due to this secondary circulation, calculated at each grid level, I, following Jeffery (2002). y ' m _ {*j(Q[<k ~ mdt] + Xjjl + l)wl+1dt} dz + [W(l + l)-w(l)]dt • [ ' where Xj are T, S, P, or N and thus the net flux of X due to the secondary circulation is W n { l ) = dtlx>V + + x ) ~ ̂ ( 0 ^ ( 0 ~ W +!) - ™(0>*j(Q] ( A 2 0 ) Appendix A. Physical Model 134 Figure A.2: Upwelling occurs over the area shown here: a fraction of the area of the ring that remains when the circle with radius r i is subtracted from the circle with radius 12- Appendix A. Physical Model Symbol Parameter Value dt timestep 15 minutes dz grid spacing 0.5 m D model depth 40 m M number of grid points 80 w x s scalar unresolved internal wave breaking shear l . O x l O - 5 m 2 s - 1 w momentum unresolved internal wave breaking shear l . O x K T 4 mV 1 vo shear instability maximum l . O x l O - 3 m2s~1 f Coriolis parameter l . l x l O - 4 r i distance from S22 to Fraser River 18.5 km distance from S3 to Fraser River 27.6 km Cyy vertical entrainment constant 3 71 vertical entrainment constant 625 mV 1 72 vertical entrainment constant 417 m 3 s - x So salinity at depth 29.53 Table A . l : Parameters used in the Physical Model. Appendix A. Physical Model 136 C l o u d Category Coefficient A Coefficient B a 0 0.4641 0.3304 0.2040 1 0.4062 0.3799 0.2344 2 0.4129 0.3420 0.2200 3 0.4263 0.3212 0.2670 4 0.4083 0.3060 0.2616 5 0.3360 0.3775 0.2787 6 0.3448 0.3128 0.2550 7 0.3232 0.3259 0.2728 8 0.2835 0.2949 0.2412 9 0.1482 0.3384 0.2381 Table A.2: Cloud model coefficients A and B (2.27) and a, the standard deviation of the trans- mission factor, T, from the regression lines. Sine of Solar Elevation Figure A.3: Cloud model fits. Atmospheric transmission, T (dots), (2.27), radiation data (Ketler, 2005) fit to T (thick line), radiation (Dobson and Smith, 1988) fit to transmission data (thin line) for cloud fraction = 1 (a) and cloud fraction = 9 (b). 137 Appendix B Biological Model IMEX scheme The biological model is of the form f i = i i j ( X ) + | { K ^ } + H / J ( X ) , 2 ) ( B . l ) where Xj are the modelled properties T, S, P, and N, Rj represent the biological processes in ((2.34) and (2.35)), Kj is the diffusive mixing term, and Wj is the upwelling term (A.20). To solve (B.I) the I M E X scheme of Ascher et al. (1995) is followed using his 7 = 1 to give semi-implicit backward differentiation X3n+l - Xjn = dt[Rj(Xn) + Dj(Xjn+1) + Wj(Xjn)} (B.2) where n represents the timestep, and Dj is the diffusion term. (B.2) defines a matrix equation 4 B n X n + i =Xf + Yl C j k (B-3) fc=i where M dimensional vectors Cjk are explicit contributions (labelled by k) from biological re- actions, Coriolis, pressure gradients, and wind stress, respectively (Jeffery, 2002). The M x M tridiagonal matrix B™ is implicit diffusion defined by equations D9a-d of Large et al. (1994). Be- cause of their inherently smaller time scale the contribution of the non-linear reaction terms Rj Appendix B. Biological Model 138 are solved using an embedded fifth-order Runge-Kutta with adaptive timestep following Jeffery (2002). Biological Tables S y m b o l Parameter Va lue Literature Value maximum growth rate depends on T [5] Ro optimal growth rate 3 . 3 x l 0 - 5 s _ 1 [1] 1 . 5 x l 0 - 5 s - 1 [6] a initial slope of uc vs. IPAR o.QSxicrVw-'s- 1 [1] Rm phytoplankton growth losses 2 . 2 x l 0 - 6 s _ 1 [2] 1 . 3 x l 0 _ 6 s _ 1 [6*] 1 . 0 6 x l 0 _ 7 s _ 1 [1*] 7 photorespiration and exudation fraction of photosynthesis 0.1 [3] K half saturation value for nitrate 0.5 fiM [4] 1 -2/zAf [6] P light inhibition parameter l . O x l O ^ W - ' m 2 [1] Table B'.l: Parameters used in the Biological Model. [1]: determined from fit to Cloern (1978); [1*]: (Cloern, 1978) value for respiration only, at 1 0 ° C ; [2]: tuned; [3]: (Jeffery, 2002); [4]: (Falkowski, 1975); [5]: determined during each time step following Eppley (1972)' [6] (Li et al., 2000), [6*] (Li et al., 2000) value for mortality only. The value R m is a combination of mortality and respiration. Addition of values for each of these processes from the literature gives a value of 2 . 3 6 x l 0 - 6 s - 1 , which compares nicely to the "tuned" value used in the model, 2 . 2 x l 0 _ 6 s _ 1 , giving us confidence that this is a realistic value. Appendix B. Biological Model 139 RUN I D Start D a t e End D a t e R l September 20 2001 November 6 2002 R2 October 8 2002 November 6 2003 R3 October 9 2003 November 11 2004 R4 October 19 2004 June 23 2005 Table B.2: Model runs used to predict the spring blooms of 2002, 2003, 2004, and 2005. First part of I D Second part of I D Third part of I D R l c R l R2 f R2 R3 w R3 R4 R4 0 half - twice (Fraser River) shift (Fraser River) one and a half (winds) 9 (cloud fraction) Table B.3: Run tests were carried out in the following way. One of Clouds (c), Wind (w), or Fraser River (f) changed and substituted with data from this table. So "RlwO" means that during R l winds were set to 0; "R2cR4" means that during R2 clouds were replaced with clouds from R4. A model test run consists of any combination of column 1, 2 and 3 of this table (except for R l with R l , etc.)


Citation Scheme:


Usage Statistics

Country Views Downloads
China 6 21
United States 3 2
Japan 2 0
City Views Downloads
Beijing 6 0
Mountain View 2 2
Tokyo 2 0
Sunnyvale 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}


Share to:


Related Items