J. Great Lakes Res. 32:852–869 Internat. Assoc. Great Lakes Res., 2006 Warmer and Drier Climates that Make Terminal Great Lakes Thomas E. Croley II1,*and C. F. Michael Lewis2 1Great Lakes Environmental Research Laboratory National Oceanic and Atmospheric Administration 2205 Commonwealth Blvd. Ann Arbor, Michigan 48105-2945 2Geological Survey of Canada (Atlantic) Natural Resources Canada Bedford Institute of Oceanography Box 1006, Dartmouth, Nova Scotia B2Y 4A2 and Graduate School of Oceanography University of Rhode Island Narragansett, Rhode Island 02882 ABSTRACT. A recent empirical model of glacial-isostatic uplift showed that the Huron and Michigan lake level fell tens of meters below the lowest possible outlet about 7,900 14C years BP when the upper Great Lakes became dependent for water supply on precipitation alone, as at present. The upper Great Lakes thus appear to have been impacted by severe dry climate that may have also affected the lower Great Lakes. While continuing paleoclimate studies are corroborating and quantifying this impacting cli- mate and other evidence of terminal lakes, the Great Lakes Environmental Research Laboratory applied their Advanced Hydrologic Prediction System, modified to use dynamic lake areas, to explore the devia- tions from present temperatures and precipitation that would force the Great Lakes to become terminal (closed), i.e., for water levels to fall below outlet sills. We modeled the present lakes with pre-develop- ment natural outlet and water flow conditions, but considered the upper and lower Great Lakes sepa- rately with no river connection, as in the early Holocene basin configuration. By using systematic shifts in precipitation, temperature, and humidity relative to the present base climate, we identified candidate climates that result in terminal lakes. The lakes would close in the order: Erie, Superior, Michigan- Huron, and Ontario for increasingly drier and warmer climates. For a temperature rise of TºC and a pre- cipitation drop of P% relative to the present base climate, conditions for complete lake closure range from 4.7T+P>51 for Erie to 3.5T+P>71 for Ontario. INDEX WORDS:Climate change, Great Lakes, hydrology, water levels, terminal lakes. INTRODUCTION Ojibway in northern Ontario and northeastern Que- bec, and drained into the Ottawa River valley and Background and Purpose thence to the North Atlantic Ocean via the St. About 9,500 radiocarbon (14C) years before pre- Lawrence River valley (Teller and Leverington sent (BP) the large upstream glacial Lake Agassiz 2004), thereby bypassing the Great Lakes basin; see in northwestern Ontario and Manitoba supplied Figure 1. melt water from the last glaciation to the upper Recent construction of an empirical exponential Great Lakes through outlets to the Lake Nipigon model of isostatic uplift for the Great Lakes region basin and thence to the Superior and Huron-Michi- following the last glaciation allowed comparison of gan basins (Teller 1985, 1987). About 8,000 14C the elevations of rebounding lake outlets with re- years BP, Lake Agassiz merged with glacial Lake constructed lake levels based on 14C-dated water level indicators such as abandoned shorelines, iso- lation basins, submerged tree stumps, and unconfor- *Corresponding author. E-mail: [email protected] mities (Lewis et al. in press a, b). The early 852 Climates Closing Great Lakes 853 FIG. 1. Paleogeography of the region north of the upper Great Lakes at 7,900 14C BP when the com- bined outflow of Lakes Agassiz and Ojibway was routed to the Ottawa and St. Lawrence river valleys (solid arrow). At this stage the water supply of the Great Lakes was no longer supplemented by inflow from upstream sources, but was supplied by precipitation alone, as at present. Prior to 8,000 14C BP, over- flow from Lake Agassiz passed into the Lake Nipigon and Lake Superior basins (open arrow). The upper Great Lakes then overflowed the North Bay outlet (open arrow) into the Mattawa River (M) and thence to the Ottawa and St Lawrence river valleys. Adapted from Figure 4p in Teller and Leverington (2004). Holocene results for the Huron, Georgian Bay, and The inferred 7,900 14C years BP low stand of the Michigan basins reveal several periods of low lake Michigan and Huron basins occurred after melt levels (lower than present) due to overflow water drainage from upstream glacial Lake Agassiz drainage through the isostatically depressed outlet began to bypass the Great Lakes basin, leaving it area near North Bay, Ontario, to the Ottawa and St. susceptible to the early dry Holocene climate (Ed- Lawrence valleys; see Figure 1. These results were wards et al. 1996). As the Great Lakes basin had anticipated on the basis of previous syntheses (Es- then entered its present hydrological regime of chman and Karrow 1985; Hansel et al. 1985; Lewis water supply by precipitation only, and as differen- and Anderson 1989; Barnett 1992; Clark et al. tial glacial-isostatic crustal uplift was accounted 1994; Colman et al. 1994a, b; Lewis et al. 1994; for, the only known process that could explain the Rea et al. 1994; Moore et al. 1994, 2000; and Lar- sub-outlet low levels is climatic reduction in water son and Schaetzl 2001). Surprisingly, however, the results showed that lake levels fell below the supply, either by enhanced evaporation or reduced Huron-Michigan basin outlet after 8,000 14C years precipitation or both. [This result for the Huron and BP, possibly for a few centuries during which low- Michigan basins is similar to the conclusions of est water levels were up to 30 m below the over- Holcombe et al. (2003) who recognized low-level flow sill at North Bay; see Figure 2. shore features beneath Lake Erie, and noting the 854 Croley and Lewis creased lake water salinity and reduced precipita- tion at the time of the closed low stands (Sarvis 2000, Blasco 2001), conditions that are consistent with reduced water supply. In Hamilton Harbour, western Lake Ontario, studies of ostracodes (De- lorme 1996), diatoms, pollen, and isotopes (Duthie et al. 1996) all reveal a low water phase and sug- gest a drier climate about 7,500 14C years BP. Simi- larly, in Mud Lake on the Keweenaw Peninsula beside Lake Superior, pollen and plant macrofossil analyses indicate an onset of a drying episode at 7,900 14C years BP with notably reduced water lev- els extending to after 7,000 14C years BP (Booth et al. 2002). As these widely-separated sites both indi- cate a lowering of lake levels and onset of a drier climatic episode, it is reasonable to envision that the inferred dry climate conditions that impacted the Huron and Michigan basins also affected other basins of the Great Lakes system. The low stand pe- riod correlates also to a relatively rapid transition in vegetative cover in the Great Lakes region from bo- real to the mixed Great Lakes-St. Lawrence Forest as recorded in pollen diagrams (McAndrews 1994, FIG. 2. Huron basin lake level (black line) Dyke et al. 2004). Complementary review and as- between 10,000 and 7,000 cal BP based on the sessment of the available paleoclimate information original elevations of lake-level indicators com- for the Great Lakes watershed, coupled with new puted by an empirical exponential model of glacial studies of proxy climatic and limnological indica- rebound for the Great Lakes basin that removes tors (pollen, isotopes, diatoms, ostracodes, and the- the effects of glacial-isostatic uplift. For clarity, camoebians) focused on the period spanning the data are removed from the plot except for the closed low stand are in progress, and will be re- Stanley unconformity constraint (S). At about ported in future publications. 7,900 14C BP (about 8,800 cal BP) the lake level, For this re-assessment, it would be helpful to ob- indicated by the Stanley unconformity (S–line tain approximate information about the amplitude indicates original elevation and circles bracket its of climatic change that might be expected to have age of 7,900 ± 300 14C BP), descended tens of caused the upper Great Lakes low stands. Accord- meters below the North Bay outlet (gray band), the ingly, we have used a hydrological model to ex- lowest possible point of overflow for the Huron plore the excursion from the present climate that basin at the time. The thickness of the grey band would force the Great Lakes into hydrologic clo- indicates the depth of water over the outlet sill at sure in terms of increased temperature and reduc- full discharge. Adapted from Lewis et al. (in press tion in precipitation. In this study, the hydrology of a, b). the Erie and Ontario basins is considered separately from that of the upper Great Lakes, as the isostati- high level of evaporative losses in the present Lake cally-depressed North Bay outlet for the upper Erie water balance, suggested that the Erie lake lakes remained at a lower elevation than the St. level may have fallen below the level of its outlet Clair River connection to Lake St. Clair and Lake sill because of enhanced evaporation sometime in Erie until much later, about 5,500 14C years BP (Es- the early to middle Holocene.] This episode of low- chman and Karrow 1985). Other attributes of the est levels appears as an extraordinarily severe im- region such as land cover, geography, and bathyme- pact of a dry climate, possibly of short duration, on try were modeled as they are at present. the upper Great Lakes hydrological system and may This low stand episode offers an opportunity, have extended to the lower lakes. once paleoclimate is better quantified, to acquire in- Preliminary study of thecamoebians and pollen in formation about the hydrological sensitivity of the the sediment sequence of Georgian Bay suggests in- Great Lakes system to high-amplitude climate Climates Closing Great Lakes 855 change. Such information would be beneficial for model studies and projections of future levels of the Great Lakes under global warming in which some climate modeling scenarios project levels below instrumentally-observed “natural variability” (Mortsch and Quinn 1996, Mortsch et al. 2000, Lofgren et al. 2002). It should be noted that the 7,900 14C years BP low stand episode occurred while the residual Laurentide ice sheet, then in the latitudes of Hudson Bay, was rapidly retreating and wasting away (Dyke et al. 2003). This was a period of rapid change in the proportions of land, ice, and water areas, with parallel changes in albedo and re- organization of atmospheric circulation (Dean et al. 2002). As a result, the inferred occurrence of closed low stand conditions in the Great Lakes basin is seen as a product of extremely unusual conditions. FIG. 3. Laurentian Great Lakes location map. It is not regarded as an analog for future conditions, but rather, as a natural experiment from which im- portant information about lake-climate sensitivity Superior, Michigan, and Huron water levels, ac- might be derived. cording to Regulation Plan 1977, under the auspices The purpose of this paper is to demonstrate that if of the International Joint Commission. a climate is extreme enough, levels on some Great Lakes Michigan and Huron are considered to be Lakes would drop sufficiently to cut off outflow, one lake hydraulically because of their connection thereby making those lakes terminal. We look at ex- through the deep Straits of Mackinac. A relatively cursions in temperature and precipitation from the small flow of Lake Michigan water is diverted into present climate to disclose those values that would the Mississippi River basin at Chicago. The water drive the Great Lakes hydrology to produce termi- flows from Lake Huron through the St. Clair River, nal lakes. This is not an attempt to simulate past hy- Lake St. Clair, and Detroit River system into Lake drology exactly, but to explore the possible Erie. The drop in water surface between Lakes magnitude of changed climates that might have pro- Michigan-Huron and Lake Erie is only about 2 m (8 duced terminal lakes about 7,900 14C years ago in ft). This results in a large backwater effect between accordance with recently acquired glacial-isostatic Lakes Erie, St. Clair, and Michigan-Huron; changes rebound evidence that Huron and Michigan basin in Lakes St. Clair and Erie levels are transmitted lake levels had descended below their overflow out- upstream to Lake Michigan-Huron. lets. From Lake Erie, the flow is through the Niagara River and Welland Diversion (used for navigation Study Area and hydropower) into Lake Ontario. There is also a The Great Lakes basin area is 770,000 km2 small diversion into the New York State Barge (300,000 mi2), about one-third of which is water Canal System which is ultimately discharged into surface; see Figure 3. The basin extends 3,200 km Lake Ontario. Lake Ontario outflows and levels are (2,000 mi) from the western edge of Lake Superior regulated in accordance with Regulation Plan to the St. Lawrence Power Project, Cornwall, On- 1958D to balance interests upstream on Lake On- tario on the St. Lawrence River. The water surface tario with those downstream on the St. Lawrence drops in a cascade over this distance some 180 m Seaway. The outflows are controlled by the Moses- (600 ft). Lake Superior is largest and deepest and Saunders Power Dam between Massena, New York has two diversions into it: the Long Lac and Ogoki. and Cornwall, Ontario. From Lake Ontario, the Lake Superior flows through the lock and compen- water flows through the St. Lawrence River to the sating works at Sault Ste. Marie and down the St. Gulf of St. Lawrence and to the Atlantic Ocean. Marys River into Lake Huron where it is joined by Lakes Superior, Michigan, Huron, and Ontario are water flowing from Lake Michigan. Lake Superior very deep (229–405 m) while Lakes Erie and St. outflows and levels are regulated to balance Lakes Clair are very shallow (6–64 m). 856 Croley and Lewis Approach We use the Great Lakes Environmental Research Laboratory’s (GLERL’s) Advanced Hydrologic Pre- diction System (AHPS), a system of hydrology, thermodynamic, and hydraulic models for the Great Lakes (Croley 2005). GLERL uses these models to make probabilistic outlooks of Great Lakes hydrol- ogy and water levels (see http://www.glerl.noaa. gov/wr/ahps/curfcst /curfcst.html), and to assess cli- mate change impacts in the Great Lakes (Croley and Luukkonen 2003, Croley 2003, Lofgren et al. 2002, Quinn and Croley 1999, Croley et al. 1998). For the purpose of this study, we adjusted the pre- sent models to simulate the Great Lakes in their pre-European-settlement natural state by removing the influences of channel control works and regula- tion plans. Also, we kept watersheds of the upper and lower lakes separate, as they were during the early Holocene, i.e., no outflow from the Huron basin to the St. Clair-Erie basin. Accordingly, we use the models here with water balances on all lakes, and with lake outflow rating curves, selected to represent “natural” or “pre-development” condi- tions. We account for lake area variations with changes in water level, but do not remove present- day diversions and consumptions as they are rela- tively insignificant for our purpose. First we consider all lakes as interdependent (as they are now, but with “natural” outlet and connect- ing channel flows) to see if simulations with histor- ical meteorology (1948–1999) produce hydrology FIG. 4. Great Lakes Meteorological Station Net- and lake levels comparable with the historical work—overland (top) and overlake (bottom). records. This allows us to assess the reasonableness of the modified models. Then, we model two sys- 1,800 stations for over-land meteorology (precipita- tems of Great Lakes independently: 1) Lakes Supe- tion and air temperature) and about 40 stations for rior, Michigan, and Huron (the upper Great Lakes), over-lake meteorology (air temperature, humidity, and 2) Lakes St. Clair, Erie, and Ontario (the lower wind speed, and cloud cover); see Figure 4. These Great Lakes) with no inflow from the upper Great data, compiled for previous studies (Croley 1990, Lakes since they drained via the North Bay outlet to Hartmann 1990, Croley 1992b, Croley et al. 1998, the Mattawa and Ottawa rivers when overflow oc- Lofgren et al. 2002, Croley 2003), provide daily curred prior to and after the low stand 7,900 14C meteorological time series over each of the 121 years ago. Next, we consider steady state hydrology riverine watersheds that drain into the Great Lakes by modeling over an extended period constructed and the seven Great Lake water surfaces. Annual by repeating the adjusted meteorological record average precipitation and air temperature are, re- until consecutive 52-year segments are identical. spectively, 80.9 cm and 2.93°C (Superior basin), We finally consider each lake as part of its parent 84.5 cm and 6.49°C (Michigan-Huron), 91.9 cm system (upper or lower system) with a water bal- and 9.19°C (Erie), and 92.3 cm and 7.41°C (On- ance on all lakes. tario). We used these historical meteorological data with our hydrology models (discussed subse- CHANGED-CLIMATE METHODOLOGY quently) to compute the “present” or “base case” The hydrology models here use daily meteoro- scenario. We then apply selected precipitation ratios logical data from 1948–1999, compiled from about and air temperature differences to the historical me- Climates Closing Great Lakes 857 teorological data and use these modified meteoro- between daily model outflows and adjusted outflow logical time series with our hydrology models to observations. Each calibration determined parame- construct changed climate scenarios. ters for infiltration, snow melt, surface runoff, per- All precipitation is adjusted by multiplying the colation, interflow, deep percolation, groundwater actual precipitation by a single precipitation ratio flow, surface storage, and evapotranspiration from and all air temperatures are adjusted by adding a all moisture storages by systematically searching single temperature difference to the actual tempera- the parameter space (with a gradient-search tech- tures. In addition, humidity is adjusted; for precipi- nique). The model agrees quite well with weekly tation ratios below unity, which are all that are and monthly outflow observations (Croley 2002, considered here, the absolute humidity is multiplied 2003). These parameters represent present-day hy- by the ratio. Thus, if precipitation is halved, then so drology and are not changed in the simulations. All is humidity. 121 model applications are used in the simulations. HYDROLOGY MODELS Evaporation GLERL’s AHPS consists of daily runoff models GLERL’s Lake Thermodynamic Model adjusts for each of the 121 watersheds, lake thermody- over-land data (original or adjusted as a changed- namic models for each of the seven water bodies, climate scenario) from the 40 over-land stations hydraulic models for the four connecting channels that are used to estimate over-water meteorology and five water body outflow points with operating for over-water or over-ice conditions based on em- plans encoded for Lakes Superior and Ontario, and pirical relationships between the two (Croley 1989, simultaneous water balances on all the lakes. It is 1992a; Croley and Assel 1994). Surface flux described in detailed overviews elsewhere (Croley processes are represented for reflection and short- 2003, 2005). wave radiation, net long-wave radiation, and advec- tion. Aerodynamic equation bulk transfer coefficients for sensible and latent heat are formu- Runoff lated with atmospheric stability effects. Energy GLERL’s Large Basin Runoff Model (LBRM) conservation accounts for heat storage; superposi- consists of moisture storages arranged as a serial tion of heat additions or losses determines tempera- and parallel cascade of “tanks” coinciding with the ture-depth profiles. Each addition is parameterized upper and lower soil zones, a groundwater zone, by age and mixes throughout the volume. Mass and and the surface channel system (Croley 2002). energy conservation account for ice formation and Water enters the snow pack, which supplies the decay. The model has been calibrated to each of the basin surface (degree-day snowmelt). Infiltration is seven lake surfaces by minimizing root mean proportional to this supply and to saturation of the square error between daily model surface tempera- upper soil zone (partial-area infiltration). Water per- tures and observations. The model enables one-di- colates from the upper to the lower soil zone and mensional modeling throughout of spatially from the lower to the groundwater zone (deep per- averaged water temperatures over the lake depth colation). Water also flows from the upper, lower, and can be used to follow thermal development and and groundwater zones into the surface channel turnovers in the lake. system, as surface runoff, interflow, and ground- water flow respectively. “Groundwater” refers to Lake Area Adjustment intra-, not inter-, watershed storage. Flows from all tanks are proportional to their amounts (linear- For each lake, precipitation p is provided as a reservoir flows). Evapotranspiration is proportional scenario-dependent boundary condition and runoff r to available water and to sensible heat (a comple- and evaporation e are estimated with the runoff and mentary concept in that evapotranspiration reduces evaporation models. They are expressed as depths available sensible heat). Mass conservation applies over the lake surface, in m, for a given time interval for the snow pack and tanks; energy conservation (day), and are based on the lake area C as coordi- applies to evapotranspiration. Complete analytical nated between the US and Canada (CCGLBHHD solutions exist. The model has been calibrated to 1977). That is, no variation of lake area is actually each of the 121 watersheds contributing to the considered in their determination in the runoff and Great Lakes by minimizing root mean square error evaporation models. However, we adjust to actual 858 Croley and Lewis lake area A by converting these depth rates into vol- cation of connecting channels through dredging or umetric flow rates, shoreline changes, use of ice control measures, and diversion of water into and out of the lakes. Any pA impacts caused by land use modification, consump- P= (1) D tive uses, and regulation of tributary rivers are viewed as reflected by changes in water supplies to the lakes and not by changes in elevation—outflow C B- A R=r (2) relationships, and were not considered in that study. B- C D We converted Southam’s relationships from their original English units and IGLD’55 water level eA E = (3) datum (CCGLBHHD 1979) to metric units and D IGLD’85 water level datum (CCGLBHHD 1995), respectively (Croley 2006). We also transformed his where P = volumetric precipitation rate in m3s–1, Lake Erie adjustment for channel project removals R = volumetric runoff rate in m3s-1, E = volumetric to one compatible with basic weir formulae and ex- evaporation rate in m3s–1, B = basin area (including pressed Ontario outflows in terms of the 1985 sill the lake), and D = number of seconds in the time in- elevation (Croley 2006). The resultant equations are terval. Note, B and C are constants for a lake while p, r, e, and A vary with time. Precipitation and Q =824.721(Z - 181.425)1.5- H , Z ‡ 181.425 (4) evaporation are directly converted by simply multi- S S S S plying the overlake rates by actual lake area. Runoff is first multiplied by the coordinated lake (cid:230) 1 1 (cid:246) 2 area (over which it was expressed) to calculate the QT =46.440Ł 2ZT + 2ZC - 166.549ł modeled runoff volume, then divided by the coordi- nated land area (to express it as the equivalent yield (cid:127) (Z - Z )12- HH , Z ‡ Z‡ 166.549 T C T T C per unit of land area), and then multiplied by the ac- (cid:230) 1 1 (cid:246) 2 tual land area to calculate the adjusted runoff vol- =46.440Ł 2ZT - 2166.549ł (5) ume. Thus “R” gets bigger as “A” gets smaller. Of course, there is some error involved with this proce- (cid:127) (Z - 166.549))12 - H , Z ‡ 166.549> Z T T T C dure since p, r, and e actually depend on actual lake area too and should have been computed from mod- els considering actual lake area and volume changes. Also, exposed land areas would not have Q =70.714(Z - 165.953)2 C C the same properties as the original basin. Consider- (cid:127)(Z - Z )-12 H , Z ‡ Z‡ 165.953 ation of the uncertainty associated with these errors C E C C E is beyond the scope of this exploratory study. =70.7144(Z - 165.953)2(cid:127)(cid:127) (6) C (cid:127) (Z - 165.953)12- H Z ‡ 165.953> Z C C C E Outflow Relations Unmanaged lake outflow depends on lake level and outflow sill elevation for lakes not affected by Q =701.504(Z - 169.938)1.-5 H , Z ‡ 169.938 (7) E E E E backwater (such as Superior, Erie, and Ontario) or on these variables as well as downstream lake level for lakes affected by backwater (such as Michigan- Q =577.187(Z - 69.622)1.5- H , Z ‡ 69.622 (8) O O O O Huron and St. Clair). (We consider the present Great Lakes here to facilitate later validation of the model.) Southam (1989) described a quantitative where Q , Q , Q , Q , and Q = outflows from S T C E O empirical relationship between water elevation and Lakes Superior, Michigan-Huron, St. Clair, Erie, outflow for each lake that represents “natural” con- and Ontario, respectively in m3s–1, Z , Z , Z , Z , S T C E ditions, prior to the introduction of societal devel- and Z = respective water elevations with respect to O opments. For the Laurentian Great Lakes the IGLD’85 water level datum in m, and H , H , S T watershed, these developments include regulation H , H , and H = respective ice retardations in C E O of outflows of Lakes Superior and Ontario, modifi- m3s-1 as shown in Table 1. Note, outflows are zero Climates Closing Great Lakes 859 TABLE 1. Great Lake outflow ice and weed when the Superior water level is below the sill of retardationa (Southam 1989). 181.425m, Q = 0 when Michigan-Huron is below T the sill of 166.549 m, Q = 0 when St. Clair is Month Superior Mich.-Huron St. Clair Erie C below the sill of 165.953 m, Q = 0 when Erie is m3s-1 m3s-1 m3s-1 m3s-1 E below the sill of 169.938 m, and Q = 0 when On- (1) (2) (3) (4) (5) O tario is below the sill of 69.622 m. January 113 1,020 425 113 Since (4)–(8) were derived from semi-empirical February 113 136 425 142 stage-fall-discharge or rating curves that were fit to March 113 651 227 85 a range of flows and elevations not necessarily April 113 170 57 142 May close to the sill, the sill elevations estimated here June 57 are in error. Sill heights on all lakes but St. Clair are July 142 well above the bottom of the lake. On Lake St. August 113 Clair, the bottom of the lake is 168.4 m (subtract September 85 maximum coordinated depth from chart datum in October 57 column 6 of Table 2); this is above the Michigan- November Huron and St. Clair sills. This corresponds to a December 113 142 channel running along the bottom of Lake St. Clair; aNo values for Ontario are given in the reference. i.e., the lake bottom is at the top of this channel and we can have flow from the Lake St. Clair basin without lake storage. Since the lake bottom is when elevation is below the “sill” elevation of the below the Erie sill of 169.938 m, we see that St. lake; the sill is the lowest elevation for which flow Clair will never be empty as long as Lake Erie is from the lake is still possible (e.g., Lake Superior’s not terminal (water line above its sill). Lake out- sill level is 181.425 m). We ignore the small eleva- flows in (4) and (5)–(8) are set to zero when nega- tion differences, introduced by the datum change, tive values would be computed (ice retardation between Michigan-Huron and St. Clair levels and would drop to equal flow rate). between St. Clair and Erie levels to keep the equa- tions physically meaningful; i.e., when Lakes Hypsometric Relations Michigan-Huron and St. Clair are at the same level The Coordinating Committee on Great Lakes (Z = Z ) or Lakes St. Clair and Erie are at the T C Basic Hydraulic and Hydrologic Data (CCGLB- same level (Z = Z ), there should be no flow be- C E HHD 1977) provided graphical relations, for each tween the respective pair of lakes (Q = 0 or Q = T C lake, between depth and volume; inspection reveals 0). However, backflow is possible from Lake Erie that simple power relations are a very good fit, to Lake St. Clair and from Lake St. Clair to Lake Michigan-Huron. This backflow is not described by V =a(M- D)b these equations (but is addressed subsequently). d (9) Note that when St. Clair water level is below the A=- V= ab(M- D)b-1 dD Michigan-Huron sill, the sill elevation is controlling in (5); likewise when Erie water level is below the where A = area of horizontal surface at depth D St. Clair sill, the sill elevation is controlling in (6). below a reference elevation, M = maximum depth, These are reasonable extensions, made here to V = lake volume beneath horizontal surface at depth allow flow computations as lake levels drop below D, and a and b are empirical parameters. By requir- those historically experienced. Note that Q = 0 ing that the coordinated values of area, C, and vol- S TABLE 2. Coordinated values of Great Lake parameters (CCGLBHHD 1977). Parameter SUP MIC HUR GEO STC ERI ONT (1) (2) (3) (4) (5) (6) (7) (8) chart datum, m 183.2 176.0 176.0 176.0 174.4 173.5 74.2 maximum depth, m 405 281 229 164 6 64 244 coordinated area, km2 82,100 57,800 40,640 18,960 1,114 25,700 18,960 coordinated volume, km3 12,100 4,920 2,761 779 3.4 484 1,640 860 Croley and Lewis ume, S, (CCGLBHHD 1977) exist at the reference C elevation (chart datum), where D = 0, for each lake, D V @ (I- QD)+ p +A R C (B - A -) e A (15) C C C C C C B - C C C C C as in Table 2, the parameters are C C a= MC D V @ (I- QD)+ p A+ R CE (B - A -) e A (16) S E E E E E E B - C E E E E (10) E E S b= Ma C D V @ (I- QD)+ p A+ R O (B - A -) e A (17) O O O O O O B - C O O O O Writing (9) in terms of elevation instead of depth, O O where D V = change in volume and the subscripts V =a(Z- Z )b B (11) refer to individual Great Lakes or extended water A=ab(Z- Z )b- 1 bodies: Superior (S), Michigan (M), Huron (H), B Georgian Bay (G), Michigan-Huron (T), St. Clair (C), Erie (E), and Ontario (O). Equation (14) con- where Z = elevation at depth D, in m, and Z = ele- B siders Lakes Michigan and Huron, including Geor- vation of lake bottom, in m, given from Table 2 by gian Bay, as one water body. Boundary conditions subtracting maximum depth from chart datum. are Water Balance I =0 (18) S The adjusted over-lake precipitation, runoff to the lake, and lake evaporation are used in a water balance, I =Q (19) T S dV =I- Q+ P+ R- E (12) I =Q , for upper Great Lakes flowing into lower dt C T =0, for upper Great Lakes flowing into the (20) where t = time, I = volumetric water body inflow Mattawa and Ottawa basins rate (outflow from the upstream lake) , and Q = vol- umetric water body outflow rate. Note that V, A, I =Q (21) and Q are not simple functions of Z. Determination E C of the proper equation to use in (4)–(8) depends on downstream water level, which in turn depends on I =Q (22) which equation is used. Likewise, backwater adjust- O E ments (described subsequently) are not reflected in (4)–(8). For these reasons, it is not possible to solve For each water body, it is necessary to compute (1)–(12) analytically. Equations (1)–(3) and (12) are the inflow as outflow from the upstream lake, the applied over time interval D to each water body lake(s) area, and the adjusted net basin supplies as based on the lakes and connecting channels part of the solution. This requires calculating lake arrangement, levels as part of the water balance. We solve (4), (5)–(8), (11) for each lake, (13)–(17), and (18)–(22) C D V @ (I- QD)+ p +A R S (B - A-) e A (13) simultaneously at each time step. Our numerical S S S S S S B - C S S S S S S procedure at each time step is: i) given p, r, e, and Z (water elevation at beginning of time step) for all D (V +V +V )@ (I- QD) 0 M H G T T lakes, ii) calculate A0 (lake area at beginning of C time step) and V (lake volume at beginning of time (cid:127) +p A +R M (B -- A -) e A 0 M M M BM - CM M M M M (14) step) for all lakes from (11) and Q0(outflow rate at beginning of time step) for all water bodies from C +p A +R H (B - A -) e A (4) and (5)–(8), iii) approximate Z (end-of-time- H H H B - C H H H H 1 H H step water elevation) as Z for all lakes, iv) calcu- 0 C late A (end-of-time-step lake area) for all lakes (cid:127) +p A +R G (B - A -) e A 1 G G G B - C G G G G from (11)and Q (end-of-time-step water body out- G G 1 Climates Closing Great Lakes 861 flow rate) for all water bodies from (4) and (5)–(8), v) approximate outflow rates and lake areas over the time increment as linear, Q = (Q + Q )/2 and 0 1 A = (A + A )/2, vi) calculate the changes in storage 0 1 for all water bodies over the time interval by using these approximate outflow rates and lake areas in (13)–(17) and (18)–(22), and vii) calculate V = V 1 0 + D V for each lake and then find Z by using V 1 1 with (11) for each water body (for Lake Michigan- Huron, interpolate for Z by using V with (11) ap- 1 1 plied to Lakes Michigan, Huron, and Georgian Bay and summed). Repeat steps iv–vii until successive values of Z for all lakes change negligibly. Repeat 1 steps i–vii for the next time step, and so forth. When solving (4), (5)–(8), (11), (13)–(17), and (18)–(22), we check and correct for backflow be- tween lakes. This could occur if water levels on Lake Erie are above those on St. Clair (and above the St. Clair sill) or those on St. Clair are above those on Lake Michigan-Huron (and above the Michigan-Huron sill). For those times when back- flow would occur between two lakes, we simply balance the lakes involved so that water levels on both are equal and the flow between them is zero. Furthermore, we consider sill heights in this adjust- ment and do not let backflow reduce a lake’s level below the upstream sill. Note that backflow does not occur when simulating the existing system with the existing climate. It also does not occur when simulating the upper lake system (Superior, Michi- gan, and Huron) with any climate since (5) is re- placed with a relation that is a function of Michigan-Huron levels only (discussed subse- FIG. 5. Net basin supply comparison for 1948- quently). Backflow corrections are only required 1999 showing reasonable agreement between when simulating the existing system or the lower water supplies based on computation of observed lake system (St. Clair, Erie, and Ontario) with lake levels and flows (historical) and supplies sim- warmer or dryer climates. The equations solution ulated from observed meteorology (modeled). converges to an insignificant difference within 2–15 iterations (the difference between water elevations els and flows). Figure 5 compares our estimates in successive iterations, summed over all lakes, is with historical NBS and shows good agreement, as less than one thousandth of a millimeter). expected since historical meteorology data are used in the simulation. Differences can be ascribed to VALIDATION water balance errors in the computation of residual To check the models and water balance approxi- NBS and to modeling errors in the computation of mations, we simulated the entire interconnected the NBS components. The biggest differences occur Great Lakes for the historical meteorological on Lake Ontario, suggesting they arise from water record. First, we compared simulated net basin sup- balance errors in computing the historical residual plies (precipitation + runoff – lake evaporation) re- NBS. sulting from the model, applied to the historical Next, we compared simulated lake levels result- meteorological record with actual initial conditions, ing from the model, applied to the historical meteo- directly to historical net basin supplies (computed rological record with actual initial conditions, as a water balance residual from historical lake lev- directly to historical levels. For this comparison, we