Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Assessing the impacts of free-stream turbines for electricity generation Atwater, Joel Fraser 2015

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

Item Metadata

Download

Media
24-ubc_2015_november_atwater_joel.pdf [ 21.42MB ]
Metadata
JSON: 24-1.0165827.json
JSON-LD: 24-1.0165827-ld.json
RDF/XML (Pretty): 24-1.0165827-rdf.xml
RDF/JSON: 24-1.0165827-rdf.json
Turtle: 24-1.0165827-turtle.txt
N-Triples: 24-1.0165827-rdf-ntriples.txt
Original Record: 24-1.0165827-source.json
Full Text
24-1.0165827-fulltext.txt
Citation
24-1.0165827.ris

Full Text

Assessing the Impacts of Free-Stream Turbines forElectricity GenerationbyJoel Fraser AtwaterB.A.Sc., The University of British Columbia, 2005M.A.Sc., The University of British Columbia, 2008A THESIS SUBMITTED IN PARTIAL FULFILLMENTOF THE REQUIREMENTS FOR THE DEGREE OFDoctor of PhilosophyinTHE FACULTY OF GRADUATE AND POSTDOCTORALSTUDIES(Civil Engineering)The University of British Columbia(Vancouver)October 2015© Joel Fraser Atwater, 2015AbstractDue to a number of factors including energy security and climate change, thereis an urgent need to transition the global energy supply to renewables. Two po-tential sources are tidal-stream and hydrokinetic power, utilizing free-stream waterturbines as generating devices.Much of the interest in tidal-stream power comes from resource assessmentsthat suggest that significant amounts of electricity could be produced from tidalcurrents flowing through straits. These assessments inventoried the kinetic energyflux and do not account for flow reduction due to turbine resistance. As such, theydo not present a realistic picture of the resource.An analytical model for flow reduction in tidal straits demonstrates that only38% of the natural fluid power is theoretically extractible. This model does notcapture the behaviour of bays, lagoons, or the open ocean. Maximum power pro-duction requires flows to be reduced to 58% of natural and if the flow is kept above95% of nominal (due to environmental regulations) less than 10% of the total poweris available.A large laboratory experiment was built to test the analytical model and theresults agree with the analytical model.iiPredicted future levelized cost of energy from tidal generation in straits is an in-terplay of reduced production due flow reduction competing with decreasing tech-nology costs. This is modelled, indicating levelized costs of energy will drop ini-tially, then rise due to flow reduction.Considering hydrokinetic power near hydropower stations, a 1D model usedSeton Canal data to simulate the installation of turbines. The results show thatthe installation of hydrokinetic turbines would decrease the output of the existingpowerhouse. Furthermore, the decrease in hydroelectric production is greater thanthe hydrokinetic production. Thus, installing hydrokinetic turbines would cause anet energy loss.In conclusion, there are three key recommendations:1. Policy makers are cautioned in embracing tidal resource assessments that arebased solely on kinetic energy flux.2. Project proponents and regulators are advised to study far-field effects of anyproposed free-stream turbine installation.3. Developers, investors, and policy makers are cautioned towards assumingthat the long-term cost of energy from tidal power will decrease.iiiPrefaceThe material included in Section 2.1 was published as Atwater and Lawrence(2011). I wrote and developed the concepts for the paper, wrote the text and madethe figures. My supervisor, Dr. Greg Lawrence, reviewed and edited the work,aided with structure, presentation and expression of the mathematical concepts.The material of Chapter 4 was submitted as a technical report to the BC Hydroand Power Authority with Joel Atwater as the sole author.ivTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiv1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21.1.1 Energy Concerns . . . . . . . . . . . . . . . . . . . . . . 21.1.2 Rate and Ease of Extraction . . . . . . . . . . . . . . . . 31.1.3 Anthropogenic Climate Change . . . . . . . . . . . . . . 41.1.4 Use of Energy and Technology . . . . . . . . . . . . . . . 71.2 Free Stream Turbines . . . . . . . . . . . . . . . . . . . . . . . . 121.2.1 Electricity from Tidal Turbines . . . . . . . . . . . . . . . 121.2.2 Electricity from Hydrokinetic Turbines . . . . . . . . . . 161.2.3 Uncertainties with Free-Stream Turbines . . . . . . . . . 191.3 Research Objectives . . . . . . . . . . . . . . . . . . . . . . . . . 201.4 Review of Literature . . . . . . . . . . . . . . . . . . . . . . . . 221.4.1 Technology Assessments . . . . . . . . . . . . . . . . . . 231.4.2 Environmental Assessments . . . . . . . . . . . . . . . . 241.4.3 General Assessment Methods . . . . . . . . . . . . . . . 252 Analytical Models for Tidal Energy Extraction . . . . . . . . . . . . 26v2.1 Regulatory, Design and Methodological Impacts on the Assess-ment of a Tidal-In-Stream Power Resource . . . . . . . . . . . . . 272.1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . 282.1.2 Model for Extraction . . . . . . . . . . . . . . . . . . . . 302.1.3 Potential Implications of Turbine/Flow Interaction . . . . 372.1.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . 432.2 Impacts of Flow Reduction and Technical Learning on LevelizedCost of Energy . . . . . . . . . . . . . . . . . . . . . . . . . . . 442.2.1 Calculation of Levelized Cost of Energy . . . . . . . . . . 452.2.2 Learning Factors . . . . . . . . . . . . . . . . . . . . . . 462.2.3 Interaction of Flow Reduction and Learning . . . . . . . . 473 Physical Modelling of Tidal Energy Resource . . . . . . . . . . . . . 553.1 Boundary Conditions and the Conceptual Model . . . . . . . . . . 563.2 Development of Physical Infrastructure . . . . . . . . . . . . . . 623.2.1 Tank Structure . . . . . . . . . . . . . . . . . . . . . . . 623.2.2 Water Supply . . . . . . . . . . . . . . . . . . . . . . . . 643.2.3 Weirs . . . . . . . . . . . . . . . . . . . . . . . . . . . . 653.2.4 Flow Resistance Element . . . . . . . . . . . . . . . . . . 693.2.5 2D Profiler . . . . . . . . . . . . . . . . . . . . . . . . . 723.2.6 Water Level Elevation Measurement . . . . . . . . . . . . 753.2.7 Bottom Resistance . . . . . . . . . . . . . . . . . . . . . 773.2.8 Electronics . . . . . . . . . . . . . . . . . . . . . . . . . 773.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 813.3.1 Calculation and Control of ∆Htide . . . . . . . . . . . . . 823.3.2 Calculation of Q . . . . . . . . . . . . . . . . . . . . . . 843.3.3 Calculation of hT and kT . . . . . . . . . . . . . . . . . . 843.3.4 Calculation of hI and kI . . . . . . . . . . . . . . . . . . . 863.3.5 Power and η . . . . . . . . . . . . . . . . . . . . . . . . 863.3.6 Operation . . . . . . . . . . . . . . . . . . . . . . . . . . 873.4 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . 883.4.1 Experimental Uncertainty . . . . . . . . . . . . . . . . . 923.5 Conclusions and Recommendations . . . . . . . . . . . . . . . . 954 Numerical Modelling of Hydrokinetic Generation Impacts near Tra-ditional Hydropower . . . . . . . . . . . . . . . . . . . . . . . . . . . 994.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1014.1.1 Existing Seton Reservoir Generation . . . . . . . . . . . . 1014.1.2 Theoretical Basis . . . . . . . . . . . . . . . . . . . . . . 1024.2 Model Development . . . . . . . . . . . . . . . . . . . . . . . . 103vi4.2.1 Selection of Model Boundaries . . . . . . . . . . . . . . . 1064.2.2 Theory and Assumptions . . . . . . . . . . . . . . . . . . 1064.2.3 Integration of Seton Canal Parameters . . . . . . . . . . . 1104.2.4 Calibration and Validation . . . . . . . . . . . . . . . . . 1124.3 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . 1134.3.1 Effects on Seton Generation . . . . . . . . . . . . . . . . 1154.3.2 Sensitivity Analysis . . . . . . . . . . . . . . . . . . . . 1174.4 Conclusions and Recommendations . . . . . . . . . . . . . . . . 1195 Conclusions and Recommendations . . . . . . . . . . . . . . . . . . 1225.1 Tidal Power . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1235.2 Hydrokinetic Power . . . . . . . . . . . . . . . . . . . . . . . . . 1255.3 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127A Literature Review of General Assessment Methodologies . . . . . . 136B Python Code for Calculation of Tidal Power LCOE . . . . . . . . . 149C Matlab Code for Control of Laboratory Experiment . . . . . . . . . 157C.1 ADV ProcessSerial.m . . . . . . . . . . . . . . . . . . . . . . . . 157C.2 ADV ProcessSerialString.m . . . . . . . . . . . . . . . . . . . . 157C.3 ADVclean.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158C.4 ADVclearLowCorrelation.m . . . . . . . . . . . . . . . . . . . . 158C.5 ADVgetdata.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 159C.6 ADVinit.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 159C.7 ADVzero.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 159C.8 calbrateUSsensors.m . . . . . . . . . . . . . . . . . . . . . . . . 159C.9 calcPower.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160C.10 CalcPowerErrorBars.m . . . . . . . . . . . . . . . . . . . . . . . 161C.11 CalcShutterK.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 162C.12 CleanUSspikes.m . . . . . . . . . . . . . . . . . . . . . . . . . . 163C.13 countdown.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 164C.14 daqinit.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 164C.15 define global.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 164C.16 DetermineFlow.m . . . . . . . . . . . . . . . . . . . . . . . . . . 165C.17 DispDepths.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166C.18 Ex7.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 167C.19 fftsmooth.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 170viiC.20 int2Arduino.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171C.21 int2string.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172C.22 mergeDAQ.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172C.23 mynanmean.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173C.24 newADVdata.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 173C.25 newDAQdata.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 173C.26 ProcessEx5.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174C.27 ProcessMultipleExperiments.m . . . . . . . . . . . . . . . . . . . 175C.28 ProcessProfile.m . . . . . . . . . . . . . . . . . . . . . . . . . . 180C.29 profileQ.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 180C.30 SetupExperiment.m . . . . . . . . . . . . . . . . . . . . . . . . . 182C.31 shutterinit.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183C.32 shuttermove.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183C.33 smooth.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183C.34 startDAQ.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184C.35 systeminit.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184C.36 travgetpos.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 185C.37 travinit.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186C.38 travmove.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186C.39 ultraClean.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186C.40 USdepth.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 188C.41 usZero.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 188C.42 weir1init.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 189C.43 weir1move.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 189C.44 weir1startserial.m . . . . . . . . . . . . . . . . . . . . . . . . . . 190C.45 weir2init.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 190C.46 weir2move.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 190C.47 weir2startserial.m . . . . . . . . . . . . . . . . . . . . . . . . . . 191D Matlab Code for Hydrokinetic Model . . . . . . . . . . . . . . . . . 192D.1 CalcPower.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 192D.2 QReduction.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 192D.3 get properties.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 193D.4 smoothClose.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 194D.5 CalculateLength.m . . . . . . . . . . . . . . . . . . . . . . . . . 194D.6 RemoveData.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 194D.7 range of TurbA.m . . . . . . . . . . . . . . . . . . . . . . . . . . 195D.8 solve flow.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 196D.9 InterpolateDepths.m . . . . . . . . . . . . . . . . . . . . . . . . . 197D.10 SensativityPlot.m . . . . . . . . . . . . . . . . . . . . . . . . . . 197viiiD.11 range of kt.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 199D.12 sub sub.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 199D.13 LoadDataFiles.m . . . . . . . . . . . . . . . . . . . . . . . . . . 200D.14 SolveKf.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 202D.15 recalc cell.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 202D.16 ModelSetup.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . 203D.17 SolveKf real.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 204D.18 recalculate channel params.m . . . . . . . . . . . . . . . . . . . 205D.19 untitled6.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 205D.20 MonteCarloTest.m . . . . . . . . . . . . . . . . . . . . . . . . . 205D.21 define global.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 206D.22 section properties.m . . . . . . . . . . . . . . . . . . . . . . . . 207D.23 PlotSensitivity.m . . . . . . . . . . . . . . . . . . . . . . . . . . 208D.24 deg2utm.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 209D.25 setup matrix.m . . . . . . . . . . . . . . . . . . . . . . . . . . . 212D.26 ProcessResultsMat.m . . . . . . . . . . . . . . . . . . . . . . . . 213D.27 depth sensitivity.m . . . . . . . . . . . . . . . . . . . . . . . . . 213D.28 smooth.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 215ixList of TablesTable 2.1 Turbine characteristics and associated flow reduction for twohypothetical tidal energy projects. . . . . . . . . . . . . . . . . 41Table 2.2 Minimum cost of energy for global resource levels. . . . . . . 52Table 3.1 Experimental measurements used to compute analytical modelvariables. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82Table 3.2 Surface elevations and ∆Htide for conducted experimental runs. 88Table 3.3 Selected data from laboratory experiments. . . . . . . . . . . . 98xList of FiguresFigure 1.1 Rate of change of radiative forcing due to emission of green-house gasses. Taken under licence from Solomon et al. (2007). 6Figure 1.2 Existing tidal power generation stations. . . . . . . . . . . . . 14Figure 1.3 Free stream turbine configurations . . . . . . . . . . . . . . . 17Figure 2.1 Extraction Ratio as a function of kT and kI . . . . . . . . . . . 33Figure 2.2 Reduction in velocity as a function of kT and kI . . . . . . . . . 34Figure 2.3 Extraction ratio as a function of reduction of flow rate fromnatural conditions. The area on the left corresponds to the re-gion where a regulatory limit on flow reduction will limit gen-eration. The area on the right represents the physical limits askT exceeds 2kI . In this region, production occurs with the flowbeing retarded less than the regulatory threshold. . . . . . . . 39Figure 2.4 Aggregation of effects that limit generation. The maximumphysical limit of 279MW for the example given is reduced byregulation and efficiency losses. . . . . . . . . . . . . . . . . 42Figure 2.5 Annual Energy Production vs Installed Turbines . . . . . . . 51Figure 2.6 Levelized Cost of Energy vs Installed Turbines. The jagged-ness of the curves is an artifact due as the model imposes aninteger number of turbines. . . . . . . . . . . . . . . . . . . . 52Figure 2.7 Levelized Cost of Energy vs Annual Energy Production . . . . 53Figure 2.8 Annual Energy Production per Turbine vs Installed Turbines . 54Figure 3.1 Constant Head Tank . . . . . . . . . . . . . . . . . . . . . . 58Figure 3.2 Conceptional model of water flume to be used. . . . . . . . . 58Figure 3.3 Rendering of the shutter plate design . . . . . . . . . . . . . . 60Figure 3.4 Arrangement of flume components. Bottom: scale of the flume.Top: detail of components. . . . . . . . . . . . . . . . . . . . 61Figure 3.5 Process flow for experiment control system. . . . . . . . . . . 63xiFigure 3.6 Two flumes in the Hydraulics Laboratory. The large yellowflume was used for this experiment. . . . . . . . . . . . . . . 64Figure 3.7 Water Supply to Apparatus . . . . . . . . . . . . . . . . . . . 66Figure 3.8 Downstream Valve and Steel Support Plate . . . . . . . . . . 67Figure 3.9 Upstream water inlet and diffuser located in upstream inflowtank . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67Figure 3.10 Downstream Weir . . . . . . . . . . . . . . . . . . . . . . . . 68Figure 3.11 Closure built for upstream outflow. Leads to in-floor trough tobe returned to the sump. . . . . . . . . . . . . . . . . . . . . 70Figure 3.12 Upstream weir . . . . . . . . . . . . . . . . . . . . . . . . . 71Figure 3.13 Shutter plate slats. . . . . . . . . . . . . . . . . . . . . . . . 72Figure 3.14 Linkages that actuate the slats of the shutter plate. . . . . . . . 73Figure 3.15 2D Traverse that supports the ADV. . . . . . . . . . . . . . . 74Figure 3.16 Ultrasonic sensor near downstream flow conditioning grid. . . 76Figure 3.17 Friction lattice on bottom of channel. Looking upstream fromdownstream weir. . . . . . . . . . . . . . . . . . . . . . . . . 78Figure 3.18 Power connections for the apparatus. . . . . . . . . . . . . . . 79Figure 3.19 Data connections for the apparatus. . . . . . . . . . . . . . . 80Figure 3.20 Cross sectional area applied to each velocity measurement forthe purposes of bulk flow calculation. . . . . . . . . . . . . . 85Figure 3.21 Change in water surface elevation along the flume as a resultof shutter closure during experimental run 3. . . . . . . . . . 88Figure 3.22 (a). Change in turbine head and flow rate as a result of shutterclosure. (b). Extracted power from flow as function of shutterangle. Data are from a Experimental Run 3. . . . . . . . . . . 90Figure 3.23 Extraction efficiency (η) versus flow rate. Each set of symbolsrepresents an experimental run. . . . . . . . . . . . . . . . . 91Figure 3.24 Beta distribution for a best guess of 4, minimum value of 3 anda maximum value of 6. . . . . . . . . . . . . . . . . . . . . . 92Figure 3.25 Variation in the shape of the velocity profile over experimentsand positions. . . . . . . . . . . . . . . . . . . . . . . . . . . 95Figure 3.26 Fraction of total fluid power extracted (η) versus turbine resis-tance relative to natural channel resistance. Error bars repre-sent upper and lower quartile of Monte Carlo uncertainty re-sults. Each set of symbols represents an experimental run. . . 96Figure 4.1 Seton Reservoir, Canal and Generation Station. Image adaptedfrom Google Earth . . . . . . . . . . . . . . . . . . . . . . . 101Figure 4.2 Effect on surface profiles after installing hydrokinetic turbines 104xiiFigure 4.3 Trapezoidal cells showing flow cross section (used to calculatevelocity) and potential head difference . . . . . . . . . . . . . 105Figure 4.4 Locations of ADCP Depth Survey Points. Image adapted fromGoogle Earth. . . . . . . . . . . . . . . . . . . . . . . . . . . 112Figure 4.5 Model and Observed Depth for 100 m3/s flow . . . . . . . . . 113Figure 4.6 Correlation of depth from model and observed depth . . . . . 114Figure 4.7 Water surface elevation at 100 m3/s flow without hydrokineticgeneration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115Figure 4.8 Water surface elevations for nominal and 64 kW of hydroki-netic generation. . . . . . . . . . . . . . . . . . . . . . . . . 116Figure 4.9 (a) Variation from natural depth as a function of distance forvarying levels of hydrokinetic generation. (b) Change in waterlevel from natural as a function of hydrokinetic generation . . 117Figure 4.10 Loss of SGS generation and gain in hydrokinetic generation. . 118Figure 4.11 Deviation from nominal depth with unadjusted and adjustedflow rate for a turbine of 200 m2 swept area . . . . . . . . . . 119Figure 4.12 Sensitivity of flow, depth and turbine drag ratio to loss of gen-eration at Seton. . . . . . . . . . . . . . . . . . . . . . . . . . 120Figure A.1 Channel flow properties with no artificial energy extraction.Recreated from Bryden et al. (2004) . . . . . . . . . . . . . . 138Figure A.2 Influence of proportional extraction on the mean flow speed.Recreated with modification from Bryden et al. (2004) . . . . 138Figure A.3 Schematic of flow through turbine or turbines, denoted by T,occupying fraction ε of cross-sectional area of channel. Recre-ated from Garrett and Cummins (2004) . . . . . . . . . . . . 140Figure A.4 Solution of Equation (A.9) for ε = {0.1,0.5,0.9} . . . . . . . 142Figure A.5 Ratio P/Pmax as function of δ1 (linear resistance). Recreatedfrom Garrett and Cummins (2004) . . . . . . . . . . . . . . . 144Figure A.6 The scaled maximum power as a function of a parameter λ ′, representing the frictional drag associated with the turbines,for various values of n where the turbine drag is assumed pro-portional to the nth power of the current speed. Recreated fromGarrett and Cummins (2005) . . . . . . . . . . . . . . . . . . 147xiiiAcknowledgmentsMy most sincere thanks to the many individuals who helped me in the process ofthis research.Thank you to my supervisor, Greg Lawrence, who first gave me the opportunityto look at renewable power, who gave me the freedom to start Crosswind Powerand HydroRun Technologies while studying, and through it always gave sound andencouraging guidance and support.Thank you to Alex Tu and Jana Hananova at BC Hydro. This research wouldnot have been possible without the contributions of BC Hydro and NSERC throughthe IPS program.Finally, thank you to Natasha, my parents, and my friends, who have alwayssupported me and all my ventures.xivChapter 1IntroductionThis dissertation is a discussion about how our society is going to meet its energyneeds in the coming decades. While myriad scholarly and mainstream articles haveattempted to address this topic as a whole, energy advancements are underpinnedby the interaction of technological, physical, economic, societal, financial, environ-mental and political forces. Due to the topic’s breadth, no single document couldcapture the nuance of all of these fields. This dissertation will focus here on themethods used in determining how much electricity may be generated using free-stream turbines in natural channels. An effort will be made to maintain connectionto the larger energy context.The motivation of this research is the recognition that a great deal of the at-tempts to inventory the tidal power resource use methodology that are fundamen-tally flawed (Atwater and Lawrence, 2011). From reading inventories such asHagerman et al. (2005), Tarbotton and Larson (2006), Polagye and Bedard (2006),policy makers, energy planners, and business leaders may obtain a very skewedunderstanding of the potential of tidal and hydrokinetic power. The following dis-1cussion will attempt to correct some prevalent misconceptions.It is important to note that the conclusions of this dissertation relate only to thebehaviour of tidal flows in straits that are driven on both ends by tidal forcing. Tidalflows into and out of bays or lagoons are driven by different boundary conditionsand are not addressed herein. Furthermore, tidal flows in the open ocean are notconsidered.1.1 Background1.1.1 Energy ConcernsThe supply of reliable, abundant, and relatively cheap energy is one of the primarydrivers for economic growth (along with labour and capital) and is thought to becritical for the functioning of western economies. Gross Domestic Product (GDP)and energy consumption are highly correlated, and there is evidence indicating thatenergy scarcity will limit economic growth. (Stern and Cleveland, 2004). Energyuse is also linked with a nation’s standard of living and overall health. Govern-ments, therefore, have established energy policies with the goal of maintaininginexpensive energy for their industries and populace.Presently, the 67% of the worlds electricity supply is derived from fossil fuel.Canada accounts for roughly 3% of global consumption, but with higher reliancenon-fossil fuel energy at 77%, compared to the global mean of 33% (World Bank,2014).There are significant risks associated with the reliance on a system that is basedheavily on the extraction of deposited oil. Non fossil-fuel sources of energy needto be developed due to the challenges surrounding the rate and ease of extraction2of the fossil-fuel deposits and anthroprogenic climate change.1.1.2 Rate and Ease of ExtractionFossil fuels are formed from the decomposition of pre-historic organic matter,which is, in essence, the storage and concentration of millions of years of incomingsolar energy. While some reformation of fossil fuel continually occurs, this rate issmall compared to the extraction rates of oil, gas and coal. As a result, the availableresource of fossil fuel is continually shrinking. At some point, these deposits willsimply be exhausted for the purposes of providing a global energy supply. (Shafieeand Topal, 2009)Another concern regarding the finite resource is the ease of extraction, particu-larly for crude oil. Oil exists in a multitude of forms, in different geological settingsand conditions. Some oil may be easily extracted simply by drilling into shallowreservoirs where the oil will exit under its own pressure, like an artesian spring.These wells tend to have the lowest associated costs. Conversely, other locationsmay be in harsh locations, requiring more intensive extraction techniques. Theseinclude enhanced recovery techniques such as injecting solvents, steam, water orcarbon dioxide into wells, or physically separating the product from surroundingmedia (such as in the tar sands). (International Energy Agency, 2014).Unsurprisingly, it is the easily accessible oil that was first extracted and is nowlargely depleted. As oil production becomes more difficult, ever increasing re-sources must be committed to extract the oil, causing a rise in the price of energy.At some point, oil production rates will reach a physical maximum and then even-tually decline; this concept is known as Peak Oil. (Boyce, 2013).In the study of petroleum geology, however, there is little consensus as to tim-3ing and sharpness of the peak of oil. Furthermore, it may be far enough in thefuture to be beyond the relevant planning horizon for most societies.Similar obstacles exist in regards to the expansion of coal extraction, thoughperhaps they are not as technically challenging to overcome as those for oil. Previ-ously, coal could be mined from concentrated seams very near the surface. Manyrecently developed mines now require the removal of large amounts of material toreach the coal deposit (Lutz et al., 2013).1.1.3 Anthropogenic Climate ChangeWhile the concept of anthropogenic climate change has been long discussed, (Ar-rhenius, 1896) it was not considered to be of great social importance until recenttimes. The first major political discussion was the UN Conference on Environmentand Development (Earth Summit) held in Rio de Janeiro in 1992, which formedthe UN Framework Convention on Climate Change.Greenhouse gasses (GHGs) are incredibly important mechanisms in regulatingthe temperature of the Earth. The planet’s mean temperature is a function of theequilibrium point of solar energy and irradiance. The amount of incoming solarradiation that is absorbed by the planet must be matched by outgoing blackbodyradiation, where the energy flux due to blackbody radiation is a proportional to thefourth-power of temperature. While there are a multitude of processes at work,greenhouse gasses are those that insulate against outgoing radiation. Higher con-centrations of greenhouse gasses will cause an increase in mean temperature suchthat the energy flux equilibrium is met.Different compounds in the atmosphere have varying greenhouse gas poten-tial due to their particular energy absorption spectra. The sources of anthropogenic4GHGs can roughly be categorized to be the result of land use change, chemical/ma-terial manufacture and energy production; energy production is, by far, the largestfraction. While there are a multitude of gasses that would cause warming if theiratmospheric concentration were increased, the primary compounds of concern dueto their anthropogenic release are carbon dioxide (CO2), methane (CH4), nitrogenoxide (N2O), sulfur hexafluoride (SF6) and several fluorocarbon gasses. Carbondioxide is by far the most prevalent of these compounds and is of particular inter-est as the majority of anthropogenic CO2 emissions are the result of the combustionof fossil fuels. (Sims et al., 2007)Since the beginning of the industrial revolution, society has caused the mass re-lease of carbon dioxide through burning coal, oil and other fuels to produce energy.Ice-core records indicate that the pre-industrial atmospheric CO2 concentration wasapproximately 280ppm; as of 2013, this concentration has risen to 397ppm (USDepartment of Commerce et al., 2014). While there is indication that CO2 concen-trations have been higher in pre-historic times, the rate of change of atmosphericconcentration is without precedence. See figure 1.1. The resultant change in globalclimate is also predicted to be without precedence and may cause catastrophic dam-age to both ecosystems and physical processes on which society depends. (Simset al., 2007)Predictions of climate impacts are made using large numerical simulations,which is a complex task that must account for the interplay between behaviouraland policy actions with physical processes. While physical process models aresomewhat deterministic, the anthropogenic inputs must be largely assumed as theyare enormously effected by policy choices, economic and population growth, andtechnological development. Multiple emissions scenarios are used as inputs to5448Palaeoclimate Chapter 6Crucifi x and Hewitt, 2005). Changes in biogeochemical cycles thus played an important role and contributed, through changes in greenhouse gas concentration, dust loading and vegetation cover, more than half of the known radiative perturbation during the LGM. Overall, the radiative perturbation for the changed greenhouse gas and aerosol concentrations and land surface was approximately –8 W m–2 for the LGM, although with signifi cant uncertainty in the estimates for the contributions of aerosol and land surface changes (Figure 6.5).Understanding of the magnitude of tropical cooling over land at the LGM has improved since the TAR with more records, as well as better dating and interpretation of the climate signal associated with snow line elevation and vegetation change. Reconstructions of terrestrial climate show strong spatial Figure 6.4. The concentrations and radiative forcing by (a) CO2, (b) CH4 and (c) nitrous oxide (N2O), and (d) the rate of change in their combined radiative forcing over the last 20 kyr reconstructed from antarctic and Greenland ice and fi rn data (symbols) and direct atmospheric measurements (red and magenta lines). The grey bars show the reconstructed ranges of natural variability for the past 650 kyr (Siegenthaler et al., 2005a; Spahni et al., 2005). Radiative forcing was computed with the simplifi ed expres-sions of Chapter 2 (Myhre et al., 1998). The rate of change in radiative forcing (black line) was computed from spline fi ts (Enting, 1987) of the concentration data (black lines in panels a to c). The width of the age distribution of the bubbles in ice varies from about 20 years for sites with a high accumulation of snow such as Law Dome, Antarctica, to about 200 years for low-accumulation sites such as Dome C, Antarctica. The Law Dome ice and fi rn data, covering the past two millennia, and recent instrumental data have been splined with a cut-off period of 40 years, with the resulting rate of change in radiative forcing shown by the inset in (d). The arrow shows the peak in the rate of change in radiative forcing after the anthropogenic signals of CO2, CH4 and N2O have been smoothed with a model describing the enclosure process of air in ice (Spahni et al., 2003) ap-plied for conditions at the low accumulation Dome C site for the last glacial transition. The CO2 data are from Etheridge et al. (1996); Monnin et al. (2001); Monnin et al. (2004); Siegenthaler et al. (2005b; South Pole); Siegenthaler et al. (2005a; Kohnen Station); and MacFarling Meure et al. (2006). The CH4 data are from Stauffer et al. (1985); Steele et al. (1992); Blunier et al. (1993); Dlugokencky et al. (1994); Blunier et al. (1995); Chappellaz et al. (1997); Monnin et al. (2001); Flückiger et al. (2002); and Ferretti et al. (2005). The N2O data are from Machida et al. (1995); Battle et al. (1996); Flückiger et al. (1999, 2002); and MacFarling Meure et al. (2006). Atmospheric data are from the National Oceanic and Atmospheric Administration’s global air sampling network, representing global average concentrations (dry air mole fraction; Steele et al., 1992; Dlugokencky et al., 1994; Tans and Conway, 2005), and from Mauna Loa, Hawaii (Keeling and Whorf, 2005). The globally averaged data are available from http://www.cmdl.noaa.gov/.Figure 1.1: Rate of change f radiative forcing ue to emissi n of greenhousegasses. Taken under licence from Solomon et al. (2007).climate models including, business as usual, and various levels of decarbonization.If the business as usual emissio s trajectory is maintained, impacts by the year2100 are predicted to include (Sims et al., 2007):1. 1-7m sea level rise2. Significant increase in storm intensity, drought, extreme temperature events3. Food shortage due to temperature changes and drought4. Ecosystem collapse5. Acidification of oceansWhile there have been a number of multinational agreements signed with re-spect to climate change, binding, global emissions targets have not yet been es-6tablished due to the social and economic consequences of reducing coal and oiluse. This may be attributed to the fact that there is no source of energy that can besubstituted at the scale, utility and low-cost of coal, natural gas or oil (US EnergyInformation Administration, 2014).1.1.4 Use of Energy and TechnologyThe problems with fossil fuels as discussed above would suggest that a prudent so-ciety would quickly find ways to move away from reliance on this energy source.While significant research is presently being conducted, there are no available tech-nologies that can substitute for fossil fuels. This is a result of numerous challengesincluding cost of energy, scalability, available resource or public acceptance.While relative proportions vary, societal energy use falls into three general cat-egories of consumption, each with its own requirements: electricity generation,transportation and direct heating. Many of the proposed strategies in reducingGHG emissions in transportation and heating involve electrification of the service,increasing the total the demand for electricity. Therefore, reductions in total GHGswill be heavily dependent on the production of non-carbon emitting electricity gen-eration.TransportationVirtually all of the energy needed for transportation is produced through the com-bustion of the petroleum derivatives of gasoline, diesel fuel, kerosene or bunkerfuel, generally referred to as liquid transportation fuels. The primary requirementof this energy source is a high energy density as the fuel reservoir needs to betransported using little space and weight.7Presently, there are two methods under development for displacing petroleumfuels for transportation: biofuels and electric vehicles.Biofuel is intended as a direct substitution for petroleum fuel to be combustedin existing or modified vehicle engines. This fuel is derived from solar energy byprocessing the cells or output of biological organisms (plants, algae or bacteria).The processing itself may require significant inputs of energy and using biofuelswill likely cause a net increase in electricity consumption. There is also concernthat large-scale increase in biofuel production could cause nutrient depletion ordisruption of the food supply.An alternative to biofuels are electric vehicles, storing energy in either a chem-ical battery or fuel cell. Electric vehicles do not represent a new energy source,only a way to store and release energy and, if deployed, will create a significantincrease in electricity demand. (Fountain, 2013)Direct HeatingA significant amount of fossil fuel energy is consumed for direct heating for spaceheating, domestic services (hot water) or industrial processes. In developed na-tions, the primary supply for this heat is natural gas or fuel oil. In developingcountries, the source may include a significant amount of biomass. (US EnergyInformation Administration, 2013b)Substitution of fossil fuel for direct heating is in the form of solar-thermal col-lection, geothermal or electrification. Offsetting heating loads will almost certainlycause an increase in electricity consumption.8Electricity GenerationStrategies to reduce GHG emissions for direct heating and transportation primarilyinvolve electrification of the process (ie. electric vehicles) or electricity intensiveprocesses (ie. heat pumps, chemical process pumps). If there is to be a reduction infossil fuel use, there will likely be a significant increase in electricity production;pure electricity demand is growing and will be further increased by shifts in trans-portation and heating. Globally, we presently produce 22,000TW-hr of electricity,the majority of which is generated using fossil fuels, primarily coal (World Bank,2014). If total fossil fuel use is to be reduced, there must be an enormous shiftin the electricity generation sources. It is important to note that if there is elec-trification of demands like heating without a generation source switch, total fossilfuel use will likely increase due to relatively low efficiencies of fuel to electricityconversion. (International Energy Agency, 2012)There are a multitude of means and a multitude of challenges in generatingelectricity without combusting fossil fuels. Presently, the most prevalent non-combusting methods of production are nuclear (11.7% of total electricity) andhydroelectric (15.6% of total electricity), however, there are significant doubtsthat these technologies can be scaled to displace the current majority incumbents.(World Bank, 2014)The generation of electricity from the nuclear fission of uranium has occurredsince 1951 when the United States developed the first electricity-producing reac-tor, (Michal, 2001). Since that time, there has been significant build-out of thetechnology with some nations such as France generating as much as 78% of theirenergy from this source (Nuclear Energy Agency, 2012). There is significant pub-9lic wariness to this technology due to significant nuclear accidents in the UnitedStates, Russia and Japan and the lack of reliable means to dispose of the nuclearwaste (US Environmental Protection Agency). As a result, the stringent permit-ting requirements for new nuclear facilities, combined with long build-times anda shortage of expertise presently limit the ability to rapidly expand the use of nu-clear technology. Significant development efforts are being expended to developfusion and thorium-fission technologies; these technologies, once developed, mayovercome some of above problems with nuclear energy (Abu-Khader, 2009).Typical hydroelectric power generation involves the construction of a large damcreating a reservoir on a river. Water from the reservoir is channelled through im-pulse or reaction turbines to convert the water’s gravitational potential energy toelectricity. Owing to the long history of hydroelectricity production, there hasbeen continual technical refinements yielding relatively low cost, dispatchable (seebelow) electricity production. Unfortunately, construction of dams and creation ofreservoirs has a number of negative environmental implications including flood-ing of land, changes in fisheries and mercury release. (Trussart et al., 2002) De-spite these challenges, a large proportion of feasible sites have already been de-veloped resulting in significant limitations on the ability to hugely expand large-hydroelectric generation.Massive reductions in fossil fuel based electricity generation will therefore belikely accomplished through the deployment of new electricity generation sources.Presently the fastest growing sources of electricity are wind turbines and photo-voltaic solar panels (Solar PV) with a number of other generation sources underdevelopment such as geothermal, wave, biomass, and the topic of this dissertation,tidal/hydrokinetic. All of these sources of generation face challenges in the current10electricity marketplace, which are examined below. With the exception of nuclearpower, all of these alternatives to fossil fuel are forms of renewable energy. (USEnergy Information Administration)Possibly the single greatest challenge to the adoption of alternative energy isthe fact there is no commercially viable or efficient means of storing large quan-tities of electricity (Nyamdash et al., 2010). As such, electricity production andconsumption must match for a given electricity distribution grid. In the presentparadigms of electricity use consumers and industry expect electricity to be avail-able on demand; this requires that electricity production must follow load. This hasbeen served by dispatchable generation, electricity production that can increase ordecrease output at a sufficient rate to match changes in load. Presently, this rapidchange in output is accomplished primarily through the use of natural gas and hy-droelectric generators (Brooks et al., 2010).Renewable energy is collected from the environment and can only be producedwhen it is available. For example, wind turbines cannot produce electricity whenthe air is still and solar panels cannot produce at night. If customers demand thatenergy always be available, ways must be found to counter the intermittency ofrenewable energy sources. This is colloquially referred to as ‘firming’ the supply.Without the development of large-scale storage technology, the leading pro-posal for firming renewable electricity supply is installing an excess and diversegeneration portfolio. The concept being that if the sun is not out, perhaps the windis blowing. Any gaps may be filled with restrained use of fossil or bio fuels. Assuch, there must be parallel development of many renewable energy technologiesbeyond wind and solar PV.Coupled with the challenge of energy storage is the challenge of maintaining11electrical grid stability. In the grid’s present form, rapid or unexpected changes ingeneration can cause a number of detrimental effects potentially causing equipmentdamage or failure of service. With generation sources that are intermittent, the sys-tem operator has no control over the generation asset’s output. For example, a gustor a lull of wind could dramatically change generation. These events exacerbatethe problem of variable consumption, therefore significant effort is presently beingspent to develop technology to condition the output of wind and solar generators.As these renewable energies collect from their local environment, there is aneed to transport the energy from the source of generation to the sites of consump-tion. This necessitates the availability of transmission infrastructure with sufficientcapacity. For example, many productive generation sites in British Columbia havenot been developed for wind due to the prohibitive cost of constructing transmis-sion (Holt and Eaton, 2008). Therefore, new renewable power generation sitesmust be sufficiently close to demand or existing transmission, or of sufficientlylarge capacity to justify the construction of transmission.1.2 Free Stream Turbines1.2.1 Electricity from Tidal TurbinesAs the moon orbits the earth, a constantly changing gravitational pull causes signif-icant volumes of water to move as waves in the world’s oceans: tides. All movingwater contains energy and it has long been realized that there is the possibility ofharnessing tides to do useful work. The earliest recorded use of tidal power sys-tems was in 787AD. (Hardisty, 2009) Energy is typically harvested from tides inone of two ways: capturing the potential energy or capturing the kinetic energy.12Tidal power has a number of favourable characteristics that address some ofthe challenges discussed above. Most significantly, this generation scheme is pre-dictable due to the precise cyclical nature of the tides. Whist it is still intermittentgeneration, electricity system operators can estimate with a great degree of cer-tainty how much energy a tidal plant will produce days, weeks or years in thefuture. This is in stark contrast to wind or solar where precise prediction cannotoccur even seconds or minutes ahead.Tidal power generation is not temporally correlated to weather or climaticevents. Tidal cycles are based only on the relative celestial positions of the earth,moon and sun. This is particularly advantageous when used in conjunction withother renewable energies as it provides additional diversity. For example, solar andwind power generation are correlated as they are both functions of local weatherwhile tidal is independent.Finally, much of the global population resides near oceans and as a result, thereis a possibility that tidal power generation can be installed relatively close to centresof population. This would reduce the potentially prohibitive cost of constructinglong transmission lines.Presently, there are three commercial tidal power plants in operation at LaRance, France (240MW), Sihwa Lake, South Korea (254MW) and Annapolis Royal,Nova Scotia (19MW). These facilities, termed tidal-barrage plants, capture the po-tential energy of the tides by impounding water in a bay at high tide. See Figure1.2a and 1.2b.At low tide in the ocean, the trapped water is released through a turbine to gen-erate electricity. While these plants have had relatively reliable operation, there aremyriad negative environmental consequences of their installation. These systems13(a) La Rance Tidal Power Station. (Public domain image.)(b) Annapolis Royal Generating Station. (Photo: Kevin McManus. Creative Com-mons Licence.).Figure 1.2: Existing tidal power generation stations.14involve the construction of a barrier across the entrance to a bay, which signifi-cantly alters fish migration and flow patterns in often sensitive estuarine environ-ments. There is increased bed and shore erosion due to the discharge jet, increasedsedimentation due to stagnant water, and increased intrusion of brackish water intothe near-shore ground water. (Gill, 2005) As a result, there appears to be littlemotivation for further development within most developed countries. A tidal bar-rage plant was proposed for the Severn Estuary between England and Wales, butapproval was denied due to environmental concerns. (BBC, 2010)As an alternative to tidal-barrage generation, much of the present technologicaldevelopment is focused on tidal-in-stream generators. These devices capture thekinetic energy inherent in tidal currents. Tidal currents are developed as large tidalwaves move about the oceans. In coastal channels, fjords, bays and straits, therising and falling of the tide can produce significant water velocities. For example,the maximum spring flood tides in the Discovery Passage, British Columbia, reach3.06 m/s (Tarbotton and Larson, 2006) and those in Pentland Firth, Scotland, reach5.7 m/s. (Shields et al., 2009)Functionally, tidal-in-stream turbines generate electricity in the same mannerwind turbines capture the kinetic energy of moving air. This class of device istermed a free-stream turbine as they interact with the flow without external struc-tures such as penstocks to induce a pressure gradient across them. In essence,free-stream turbines rely on the inertia of the working fluid to force flow throughthe turbine. A great number of design configurations for tidal-in-stream turbineshave been proposed. These include:• Open blade horizontal axis (much like conventional wind turbines)15• Ducted horizontal axis• Vertical axis (typically a Davis, Darius or Gorlov turbine)• Cross flow• Counter rotating• Oscillating hydrofoil• Biomimetic finsThe devices furthest along the development pathway appear to be those of the un-ducted and ducted horizontal axis configurations, as shown in Figures 1.3a and1.3b. (Khan et al., 2009)It has been suggested that these turbines will have a less acute environmentalimpact than tidal barrage plants due to a smaller footprint owing to the absence oflarge impoundment basins. (Gill, 2005) Presently, there have been only a handfulof relatively small tidal-in-stream devices in operation. As a result, the predictionof relatively benign environmental impact has not yet been tested at larger scalesof generation.1.2.2 Electricity from Hydrokinetic TurbinesHydrokinetic turbines are very similar in configuration to typical tidal turbines. Thedifference between the two classes of devices is primarily not with form but withthe environment in which they operate. Hydrokinetic turbines typically (and in thisdissertation) refer to free-stream turbines used in rivers or canals and are contrastedto hydroelectric, large-hydro, or run-of-river hydro (termed conventional hydro) bythe absence of pressurized penstocks.16(a) Characterized unducted horizontal axis turbine.(b) Characterized ducted horizontal axisturbine.Figure 1.3: Free stream turbine configurations17Similar to the difference in environmental impact between tidal-barrage andtidal-in-stream generations, it is thought that hydrokinetic generators will not dam-age the local environment or fish populations as critically as large-hydro or run-of-river hydro, (Amaral et al., 2011).While both tidal-in-stream and hydrokinetic turbines operate using identicalmeans to extract kinetic energy from moving water, the forces causing that waterto move are dramatically different. Flow in rivers is usually the result of rainfall,groundwater discharge, or snow/glacier melt while canal flows may also be theresult of water pumping. As a result of its source, hydrokinetic generation doesnot experience the same level of intermittency as tidal power. Excluding opera-tional changes, power generation is a function of river/canal velocity. Dependingon the water source and the hydrology of the watershed, there tends to be seasonalvariation in basal flow with spikes and lulls due to rain and weather events. If ahydrokinetic turbine is running at full capacity due to seasonal basal flow, electric-ity system operators can have some confidence that there will not be large changesin the output of the generator. Furthermore, changes in river flows tend to occurin the time-scale of hours, rather than seconds of wind or solar power. The longertime scale of hydro results in greater ease of electric grid control.While a large fraction of the human population resides near oceans, the ma-jority reside near rivers due to their need to access fresh water. (Kummu et al.,2011) While it is not always feasible to construct conventional hydro near cities,hydrokinetic has the potential to generate electricity very near and within centresof population.181.2.3 Uncertainties with Free-Stream TurbinesLike all new technologies under development, there are a number of uncertaintiesthat exist. These uncertainties can affect the feasibility of the technology to be de-ployed for widespread use. Such risks, if categorized in a way common to businessanalysis, include the following:Technological Risk:• Can the devices survive in the marine environment? (Examples: corrosion,debris, bio-fouling)• What is the operating lifetime of the device?• Can the production be scaled up?• What is the power quality produced by the devices?• Can the devices be easily installed and connected?• Could these devices harm the environment or marine life?Market Risk:• Will the cost of energy produced be below grid price?• Does the technology rely on incentives? If so, are these incentives perma-nent?• Is there sufficient resource to build-out the technology?• Will system operators purchase electricity produced from tidal/hydrokineticgeneration?19• Can productive generation sites be identified?• Will the predictable generation of tidal power give a premium price for elec-tricity?• Can other technological advancements (such as storage) eliminate the needfor diversity in renewable energy?Regulatory Risk:• Will these devices be permitted in navigable waters?• What are the acceptable impacts on fish or wildlife?• What are the acceptable impacts on natural water movements?• What bonding/clean-up requirements will be imposed on operators?• Will the devices impede other users of the water body/water way?• What certification standards will be applied to the devices?Clearly there is a multitude of uncertainties to be addressed. Some of theserisks are minor or may be mitigated through sound engineering. Others are randomor can be mitigated, but cannot be fully determined until policy makers or marketsaddress the issues. There are some risks, however, that exist due to lack of researchand these, will be the primary topic of investigation in this dissertation.1.3 Research ObjectivesThe risks identified in the Section 1.2.3 apply to a number of groups includingtidal/hydrokinetic technology developers, utilities and governments. The objective20of this research is therefore to identify and mitigate some of the risk by reducinguncertainty and doing so in a way that is meaningful to the above groups.From the potential risks to the tidal and hydrokinetic industry, the followingtwo, broad questions will be questions addressed:• Can we determine if there is sufficient resource to build-out the technologyas an industry?• How can the impacts of generation on natural flows be understood?The two research questions addressed in this dissertation are independent ofthe actions of the entities involved in the tidal/hydrokinetic sector. For example,technology developers may design their turbines to withstand the impact of de-bris; regulators dictate the acceptable levels of fish mortality through turbines; orutilities may choose to, or not, pay premiums for predictable power.In contrast, the available resource for tidal generation, for example, may bedetermined by inspecting the flow physics. However, the goal is not to producemodels that predict flow precisely and answer these questions exactly. The goal isto resolve these questions well enough to allow industry, government and utilitiesto reduce the uncertainty and risk. Importantly, these questions must be answeredin a manner that is accessible and appropriate to these parties.While a comprehensive global or regional inventory of tidal and hydrokineticsites is beyond the scope of this dissertation, methods will be developed for entitiesto assemble this inventory for their region of interest. Examples of the applicationsinclude allowing the following:• businesses to determine their total addressable market;21• utilities and system operators to plan future balances of intermittent, base-load dispatchable generation assets; and• governments to determine the resource so they may:– set levels of incentive to achieve a particular energy mix,– allocate appropriate funding to feed-in-tariff programmes,– plan their energy portfolio and future greenhouse gas projections.Given the existing literature on the topic, the following hypotheses are thereforepresented:• Tidal resource can be limited by regulated or physical flow reduction;• Hydrokinetic generation can cause significant impacts on existing generationassets in a waterway; and• The results of the analytical models that have been developed will be sup-ported by physical experiments.1.4 Review of LiteratureIn keeping with the general interest in the field of tidal energy, there are a numberof academic articles addressing the topic. Upon review of many of these arti-cles, they may be grouped into two broad categories, technology assessments andenvironmental assessments. Technology assessments detail the performance andbehaviour of turbine and its immediate surrounding area, while environmental as-sessments are concerned with the interaction of the turbine and the larger system.221.4.1 Technology AssessmentsMany technology assessment publications are either characterizing the performanceof a particular turbine to a variety of conditions, or characterizing turbine to tur-bine interactions. An example of turbine characterization can be found in Kim et al.(2012). These papers are specific to the turbine geometry under consideration andconsider the flow the turbine experiences to simply be a boundary condition of theexperiment or computational model.The other broad class of technology assessment articles addresses turbine toturbine interaction. If tidal turbines are to be installed in clusters with multipleunits, there may be interactions between the devices, largely due to wake. A num-ber of articles such as Myers and Bahaj (2010) or Vennell (2012) have studieda variety of turbine arrangements and spacings using experimental and computa-tional methods. Once again, these studies do not consider the feedback of turbineresistance on the tidal flow, considering it only an input.Finally, there are a small number of review articles such as Khan et al. (2009)or Rourke et al. (2010) which outline the variety of styles and attributes of dif-ferent turbine types. There is still considerable development being undertaken intidal turbine technology, and the industry has not converged on a particular designconcept.As was discussed above, the research goals are to determine the magnitude ofthe tidal energy resource and its interaction with flows. While there must be anunderstanding of the general performance characteristics of turbines, the specificparameters may be reduced to a range of performance and thrust coefficients ofturbines, per general theory regarding turbo machinery.231.4.2 Environmental AssessmentsThe articles addressing environmental assessments are comprised of four generalgroupings: impacts of turbines on biota, inventories of kinetic energy flux for aregion, modelling of tidal energy extraction for a region, and tidal energy resourcemethodology. Only the last of this this set of is relevant to this dissertation.• Biota Impacts: The articles discussing the interaction of turbines and biota,such as Gill (2005) and Clark (2006) discuss the range of impacts includingfish and seabird mortality resulting from turbine blade strikes. Due to thesmall number of long-term tidal power installations, the data is somewhatlimited. These studies may be valuable to regulators who are tasked withpermitting tidal power developments. The behaviour of marine ecosystemsis largely superimposed on the physical system however. Changes in thephysical system can impact the ecosystem, but changes in the biota are un-likely to have an effect on flows as energetic as those found in tidal channels.• Kinetic Energy Flux Inventories: There are many articles such as Tarbottonand Larson (2006), Hagerman et al. (2005), Polagye and Bedard (2006) thatpurport to be assessments of tidal energy resource, but their methodology issimply an inventory of kinetic energy fluxes in a region. It has been shown byAtwater and Lawrence (2011) (See Section 2.1), and Garrett and Cummins(2005) among others that kinetic energy flux is not related to the tidal energyresource. It may be surmised that the large number of these publicationsare the result of the endorsement of the method by the utility group EPRI(Hagerman et al., 2005).• Tidal Energy Assessments for a Region: Recognizing that the installation24of turbines will have an impact on tidal flows, studies have assessed theresource in a region using computational methods, such as Walters et al.(2013), Draper et al. (2014) or Ramos et al. (2013). Typically, a conventionalcomputational code is used when the resistance caused by the turbines issimulated as artificially high bottom friction. This type of simulation is alsovery useful in determining optimum placement of turbines within a channel.The simulation is dependent on the specific parameters of the location, sowhile useful for planning of that location, the results of these articles do notaid in the discussion of the global resource or tidal energy implications.• Resource Assessment Methodologies: There has been relatively little workdeveloping general assessment methodologies for tidal energy resources asidefrom computational simulation of a region. While useful, general conclu-sions cannot be drawn from these computational publications. The authorof this dissertation performed a review of the available literature in Atwater(2008). Since that time, no additional general resource methodology papersof note have been found. The relevant sections from Atwater (2008) havebeen included in Appendix A for reference.1.4.3 General Assessment MethodsThere are a number of publications that discuss methodology to assess tidal powerresource such as Bryden et al. (2004), Garrett and Cummins (2004) and Garrettand Cummins (2005). These articles were reviewed in the author’s master’s thesisand are included in Appendix A.25Chapter 2Analytical Models for TidalEnergy ExtractionIt is the inevitable consequence that if energy is produced from a natural system,that system will be altered. The term ”produced” is of course a misnomer, yet inthe parlance of renewable energy ”produce” is used interchangeably with the moreaccurate terms ”extract” and ”convert”. The first law of thermodynamics excludesthe possibility of any system from producing energy; it may only be converted fromone form to another. Thus, if the output of natural system under consideration hasa new output of electricity that wasn’t present in the natural state, the other energystores or flows of the system must have changed.The first question to determine is whether a meaningful amount of energy canbe extracted and converted without altering the system to such an extent that itsability to continue to deliver energy is compromised. The term meaningful is sub-jective and dependent on context. For the purposes of this dissertation, a resource26will be considered meaningful if it could be developed to a point where it providesa substantial portion of a regional energy supply.We are not attempting to ascertain if developing the resource is prudent. In-stead, the goal is to determine if it is even feasible within the general constraints andexpectations of reasonable technological development, energy markets and regula-tory environments.There have been a number of studies that suggest a meaningful amount of en-ergy could be extracted from tidal flows using free-stream turbines. These studies,however, have almost universally used a methodology, kinetic energy flux, that isfundamentally invalid for this purpose. The following section 2.1 discusses theshort-comings of the that methodology, it presents a new model and discusses theregulatory and economic impacts. Section 2.1 was published as a article in theJournal of Energy Policy (Atwater and Lawrence, 2011).Using the analytical model presented in section 2.1, relationships are built be-tween total electricity generated from ocean tides and combined with a technologyadvancement model, presented in section 2.2. The change in levelized cost of en-ergy from tidal power is examined as the technology progresses through its rollout.2.1 Regulatory, Design and Methodological Impacts onthe Assessment of a Tidal-In-Stream Power ResourceTidal-in-Stream energy has been heralded by many as a significant potential sourcefor clean power, a scheme where kinetic energy is extracted from tidal currents.A number of estimates have suggested that tidal power may become a sizeablefraction of overall electricity generation, however these estimates have been largelybased on a resource assessment methodology that dramatically oversimplifies the27physical phenomenon at play. We develop a model that considers the effect ofenergy extraction on the bulk flow of tides in a channel, showing that tidal energyinventories that assess solely kinetic energy flux may represent both an order ofmagnitude overestimation of the resource and a significant oversimplification ofregulatory impacts. Tidal flows associated with bays, or tidal flows in the openocean unconstrained by land are not considered.The interplay between the characteristics of a flow and the regulatory and eco-nomic issues will likely limit tidal power generation to levels significantly belowthe physical maximums. Permitted flow reduction, turbine design and staging ofdevelopment all have significant and predictable impacts on the extractible re-source. Energy planners must therefore understand these relationships in orderto appropriately assess the magnitude of generation that can be realistically be pro-duced from tidal energy.2.1.1 IntroductionAs a result of public pressure and growing legislation aimed at mitigating the ef-fects of global climate change, utilities around the world are looking for new meansof generating electricity without the associated release of greenhouse gas. A myr-iad of technologies has been proposed to renewably generate electricity, the mostvisible being wind and solar, however these technologies have problems in both in-termittency and unpredictability. Given that there are presently no feasible meansfor large scale energy storage, it is felt that diversification of generation beyondwind and solar is necessary to ensure a secure energy supply (Goldemberg et al.,2000).One such alternative generation scheme is tidal power. The kinetic and poten-28tial energy of the tides has long been utilized: two major tidal-barrage power plantsare presently in operation: Annapolis Royal, in Nova Scotia (built in 1984, generat-ing 18MW) and La Rance, France (built in 1966, generating 240MW) (Hammons,1993). These facilities have operated for several decades, however some signifi-cant local environmental impacts such as sea-bed erosion, fish kills, groundwatersalination, and disruption of the estuarine system have been observed (Tidmarsh,1983). As a result of the environmental impact, there has been little further de-velopment using this style of plant. The proposed barrage in the Severns Estuaryin Wales was cancelled for this reason (Berkowitz, 2010), however a handful offacilities in South Korea are under development (Bae et al., 2010).A newer tidal scheme, termed Tidal-in-Stream Energy Conversion, involvesplacing free-stream turbines in a channel where tidal fluxes create significant cur-rents. While the environmental impact is somewhat difficult to ascertain due to thelimited development thus far, it is thought that tidal-in-stream devices will be sig-nificantly more benign for marine life when compared to tidal-barrage plants (Gill,2005).Even as companies are developing turbine technology, there is significant un-certainty as to the limits of the tidal-in-stream resource. A number of assessmentshave been completed that simply take an inventory of the kinetic energy flux over achannel cross-section (Triton Consultants Ltd., 2002; Tarbotton and Larson, 2006;Polagye and Bedard, 2006; Previsic et al., 2006; Hagerman and Bedard, 2006;Hagerman et al., 2006). Implicit in the use of this methodology is that the bulkflow in the channel is unaffected by the introduction of turbines. While this as-sumption may be reasonable for very small levels of extraction, it is invalid whenassessing the limits of the total resource.29Addressing the problem of kinetic energy flux inventories being fundamentallyinappropriate, Garrett and Cummins (2005) developed a more robust methodol-ogy for determining the extractible power. This involves the integration of theone-dimensional momentum equation along the length of the channel. While thisapproach provides a defensible estimation of the available power, it may offer lim-ited utility for those wishing to extend the analysis in the context of a resourceassessment. In this paper, we develop a methodology in parallel with Garrett andCummins (2005) which can be applied with limited environmental data, but is ex-tensible and allows for investigations of the practical implications of tidal energydevelopment.2.1.2 Model for ExtractionHead Loss Along a ChannelConsider a tidal channel, each end open to the ocean at distant points. As a tidalwave travels past this channel, a head difference is established between the twoends which drives a flow. The proportionately larger size of the ocean relative tothe channel provides that the flow will not have an effect on the ocean water levels.Furthermore, the period of the long wave traveling through the channel is muchshorter than the 12.5 hour forcing period of the tide and the system may thereforebe considered quasi-steady.Analysis of tidal flows associated with bays, lagoons, and the open ocean areexcluded. Sites such as Skookumchuck Narrows present high tidal velocities, butdo not represent the majority suitable the sites as identified in Tarbotton and Larson(2006), or globally significant sites such as Pentland Firth. As such, these sites are30neglected in the analysis of this dissertation.The flow in the channel is driven by a hydrostatic pressure differential andopposed by a combination of natural bottom friction and the resistance providedby any installed turbines. This relationship may be written in terms of head lossesas∆Htide = hI +hT +u22g(2.1)where ∆Htide is the difference in tide level between the ends of the channel, hI andhT are the head losses associated with bottom friction and the turbines respectively,g is gravitational acceleration and u is the mean fluid velocity in the channel. Theinclusion of the velocity head term ( u22g ) on the right-hand-side of Equation (2.1)accounts for the flow leaving the channel downstream as a jet and dissipating inthe ocean. Given that this term accounts for a frictional process, it is unnecessaryfor it to remain independent and may be included within hI .Head loss due to friction is considered in the conventional form, as being pro-portional to the square of the fluid velocity allowing kT and kI to represent headloss coefficients associated with hT and hI respectively. Equation (2.1) can now beexpressed as∆Htide = (kT + kI)u2 (2.2)Bulk velocity in Equation (2.2) is then solved for, givingu =√∆HtidekT + kI(2.3)Thus, the velocity can be determined by considering only the tidal head difference(∆Htide) and the resistance coefficients of the turbine (kT ) and bottom friction (kI).31Generated PowerThe power extracted by a turbine may be considered as the rate the fluid can dowork on that turbine as given byPextracted = ρghT Q (2.4)where ρ is fluid density and Q is flow rate. While Q in Equation (2.4) relatesonly to the flow incident to the turbine, Garrett and Cummins (2004) determinedthat the maximum extraction occurs when the turbines occupy the entire channelcross-section, therefore for the purposes of assessing resource limits, Q in Equation(2.4) may be considered to be the product of the bulk velocity and channel cross-sectional area.This leaves the extracted power, a function of the tidal head difference, naturalresistance and resistance of installed turbines asPextracted = ρgAkT(∆HtidekI + kT)3/2(2.5)Extraction RatioIt is useful now to define a new term, the extraction ratio η , as the relative mag-nitude of the extracted power in relation to the natural fluid power of the channelsystem, Pnatural such thatη =Pextractedρg∆HtideQnatural(2.6)where Qnatural is the flow rate through the channel before the installation of tur-bines. Combining Equations (2.3), (2.5) and (2.6), the extraction ratio can be ex-32pressed asη =kT√kI(kT + kI)3/2(2.7)This relationship demonstrates that a local maximum exists with respect to the re-sistance of the installed turbines as shown in Figure 2.1. Differentiating Equation(2.7) with respect to kT and setting the result to zero gives the conditions of maxi-mum extraction ratio asηmax =233/2= 0.38 kT = 2kI (2.8)indicating that the physical limit to extraction is no more than 38% of the fluidpower of the channel. Note, that this limit is independent of the Betz limit thatrelates to the energy extracted by the swept area of an un-ducted turbine(?).0 5 10 15 2000.10.20.30.4Relative Turbine Resistance (kT/kI)Extraction Ratio (η)Student Version of MATLABFigure 2.1: Extraction Ratio as a function of kT and kI33Response of Flow RateFor any non-zero value of kT , there will be a reduction in the tidal flow, the relativemagnitude of which can be found from Equation (2.3) asQQ0=√kI√kT + kI(2.9)The rapid reduction in flow rate is shown in Figure 2.2. By applying the conditionsfrom Equation (2.8), the flow rate at maximum production will beQQ0=1√3= 58% (2.10)relative to the natural, undisturbed flow.0 5 10 15 20 25 3000.20.40.60.81kT/kIQ/Q 0Figure 2.2: Reduction in velocity as a function of kT and kI .Temporal Variation of Tidal HeadAnalysis to this point has been conducted considering a constant tidal head differ-ence but may be applied to a time-varying tide in two ways: (i) Evaluating a timeseries as discrete values for ∆Htide and solving Equation (2.5) for each time-step34or (ii) Approximating the tidal forcing as a sinusoid. If the time varying tide isexpressed as∆Htide(t) = ∆Htide︸ ︷︷ ︸maxsin(ωt) (2.11)where ω is the frequency of the tidal forcing, a capacity factor may be given asP¯Pmax=∫ 2piω0 (|sin(ωt)|)3/2 dt2pi/ω= 0.556 (2.12)This factor can be multiplied by η , giving the average power over a tidal cycle asP = 0.221ρg∆Htide︸ ︷︷ ︸maxQ0︸︷︷︸max(2.13)This solution agrees with Garrett and Cummins (2005) for the case of friction dom-inated quadratic drag and gives a order-of-magnitude assessment of the resourcewith only limited environmental data required.Example Idealized ChannelWe now apply the previous methodology to investigate the implications of flowreduction. Consider a very significant tidal channel that is 1.5km wide by 40mdeep with a maximum velocity of 5m/s and a tidal range of 4m where the tideis delayed 25 minutes over the length of the channel. Note that while this is anidealized channel, the geometric sizing loosely approximates Discovery Passage,British Columbia, Canada.The height of the tide is sinusoidal such that the water surface elevation at eachend of the channel is given by α1 sin(ωt) and α2 sin(ωt−φ) where α1, α2 = 2m.Harmonic addition theorem provides the difference in these sinusoids will result in35a third sinusoid, the amplitude of which represents ∆Htide, as given by∆Htide =√α21 +α22 −2α1α2 cos(φ) (2.14)This gives the maximum ∆Htide = 0.42mThe maximum flow rate is given simply by the product of the maximum ve-locity and channel cross-section, such that Qmax = (5m/s)(40m · 1500m) = 300 ·103m3/s.From equation (2.13), the maximum extraction isP = (0.221)(1025kg/m3)(9.81m/s2)(0.42m)(300 ·103m3/s) = 279MW.Let us now compare this prediction of tidal power using kinetic energy flux, asgiven by:PKE =12ρAu3 (2.15)Similar to Equation (2.16), average kinetic energy flux can be found usingP¯Pmax=∫ 2piω0 (|sin(ωt)|)3 dt2pi/ω= 0.424 (2.16)This yields the mean kinetic energy flux asPKE = (0.424)(0.5)(1025kg/m3)(40m ·1500m)(5m/s)3 = 1630MWComparing to the results that consider flow reduction, kinetic energy flux over-estimates the resource in this example by nearly six-fold.36Discrepancy of Kinetic Energy FluxThe ratio between extractible power and kinetic energy flux is a function of the umaxand ∆Htide. From equations (2.13) and (2.15), the general relationship betweenextractable power and kinetic energy flux is found to be:PextractedPKE=g∆Htide0.96u2max(2.17)High velocities are required for tidal power production to be technically andeconomically feasible, thus the locations of interest would be those with strongtidal currents (Tarbotton and Larson, 2006). As such, the appearance of the u2maxterm in Equation (2.17) raises significant practical concern as ‘ideal’ generationlocations (those with a large u2max) may tend to be biased towards kinetic energyflux overestimating the resource. Some insight may be gained by investigating ifthe large velocities are caused by a large head difference driving the flow (large∆Htide) or a constriction temporarily accelerating the flow without a correspondingpotential head.2.1.3 Potential Implications of Turbine/Flow InteractionOur analysis above demonstrates the physical limits on tidal power generation,but we have not yet considered practical limits due to environmental regulations,turbine design or economics. It will likely be these factors that constrain tidalpower generation.37Factor I - Environmental RegulationIn order to produce the maximum extractible power described as in Equation (2.8),the tidal flow must be reduced to 58% of natural. Given the importance of tidalcurrents on estuarine dynamics, effluent flushing and marine life, it may not beprudent for western environmental ministries to allow such a dramatic reduction.To estimate the impact of regulation, equations (2.7) and (2.9) are combined togive the extraction ratio as a function of permitted flow reduction asη =QQ0−(QQ0)3(2.18)The relationship between flow rate reduction and extraction ratio is shown in Figure2.3. If the regulatory limit requires the flow to remain at some point above 58% ofnatural, this regulation will act as the limit of generation. However, if the regulatorylimit is below 58%, it will have no effect as the physical limits will constrain powerproduction.While it is difficult to speculate on exactly what a particular government maypermit, reductions on the order of 5% are discussed in the literature as a potentialbenchmark (Hagerman et al., 2005). This allowable reduction provides that lessthan 10% of the fluid power may be extracted (contrasted to the physical maximumof 38%).Factor II -Equipment Performance and ArrangementAll discussion of energy extraction previous to this point has focused on loss ofenergy from the fluid, not the generation of usable energy. There are two factorsthat cause the energy extracted to differ from the usable energy produced, drag on380 20 40 60 80 10000.050.10.150.20.250.30.350.4% Permitted Flow Reduction (1ïQ/Q0)Maximum Extraction Ratio (d)Maximum d.Generation limited by physical potential.Generation limited by regulation.d if flow were reduced to full permitted extent.Figure 2.3: Extraction ratio as a function of reduction of flow rate from nat-ural conditions. The area on the left corresponds to the region wherea regulatory limit on flow reduction will limit generation. The area onthe right represents the physical limits as kT exceeds 2kI . In this region,production occurs with the flow being retarded less than the regulatorythreshold.auxiliary structures and the efficiency of the turbine, generator and transmissionsystems.The turbine resistance, kT , represents the total drag applied by the sum of theturbine blades, hub, ducting, support structure, anchors, etc. Drag on the structurewill act to retard the flow but not generate useable power; similarly an inefficientturbine will dissipate energy from the flow (through fluid friction on the blades)without usable power production.This introduces a potential conflict in priorities between turbine developers andregulatory agencies in relation to the design of tidal energy devices. Consider aug-39mentations such as ducting or shrouding around a turbine rotor. Adding ductingcan increase the power output of a device by allowing the turbine to intercept across-section of flow larger than its swept area. While this may improve the eco-nomic performance of the device, the auxiliary drag associated with the ductingwould have the effect of diminishing the overall extractible energy. The turbine de-veloper would favour the lower costs of generation, while regulators would likelybe concerned with maximizing the overall use of a finite tidal power resource.When issuing permits, governments must therefore recognize the difference be-tween licensing the level of electricity/energy production versus the level of energyextraction.Factor III -Reduction of Power DensityIt is important to note that incremental addition of resistance that reduces the tidalflux has effects beyond new devices. All turbines in the channel intersect the sameflow and therefore experience the same reductions in that flow. As the ability toextract energy is dependent only on the turbine’s characteristics and the cube of thefluid’s incident velocity, the performance of previously installed turbines will dropas new units are installed.This scenario raises potential disputes between developers and between devel-opers and policy makers. As an hypothetical example, consider Developer A andDeveloper B with project characteristics as shown in Table 2.1. Developer A in-stalls turbines with performance specifications to match the predicted flows. Atsome later point, Developer B installs its project which causes the total volumetricflux to drop 4.5% (to 93% of natural). Due to the cubic term in turbine perfor-mance, Developer B’s installation would cause Developer A’s energy production40Table 2.1: Turbine characteristics and associated flow reduction for two hy-pothetical tidal energy projects.Initial Installation Subsequent InstallationDeveloper A kt 0.05 0.05Developer B kt - 0.10Total kt 0.05 0.15Q/Q0 97.6% 93.3%Fluid Power Density 100% 87.2%to fall by 13%, coupled with the associated economic consequences. This situationis exacerbated as turbines are most efficient when they are operated at the designspeed and in changing the fluid velocity, the existing turbines will have sub-optimaloutput.This effect grows as the limits of the resource extraction are reached whenper-generator output falls to 19% of initial conditions (potentially producing afive-fold increase in costs). Royal Academy of Engineering (2004) predicts thatwhen tidal energy generation technology is mature but without considering flowreductions, tidal-in-stream generation costs will be approximately 0.057GBP/kWh(US$0.116/kWh). These rates are only marginally economically viable and a five-fold increase in cost would likely make tidal power production prohibitively ex-pensive.Aggregation of Secondary EffectsConsider again the idealized channel as was presented above, now with the lim-itation of 5% flow reduction, a typical turbine efficiency of 30%, a structure thatcauses drag equivalent to the turbine itself and generator and transmission efficien-cies of 90% and 80% respectively. Beginning with the physical limits of extraction(279MW), each of the above factors compound to reduce the amount of electricity41that may be produced. From Equations (2.18) and (2.6), energy generation wouldbe limited to 7.2MW representing 2.6% of the physical maximum and 0.44% ofthe kinetic energy flux as shown in Figure 2.4. In this example, the kinetic energyflux (as it would be listed in a typical resource assessment such as Hagerman et al.(2006)) provides a value that is more than two-hundred times larger than the powerthat is realistically extractable.050100150200250300TotalExtractible PowerFlow ReductionLimitedTurbine Efficiencyand DragGeneratorEfficiencyTransformer andTransmissionEfficiencyElectricityGeneratedPower (MW)279 (211)(2.81)(1.00)(57.0)7.24Figure 2.4: Aggregation of effects that limit generation. The maximum phys-ical limit of 279MW for the example given is reduced by regulation andefficiency losses.422.1.4 ConclusionsWhile Tidal-in-Stream power may appear to be a potentially viable source forlarge-scale energy generation, the assessment methodology previously used to in-ventory the resource is fundamentally inappropriate. Bulk kinetic energy flux doesnot account for the effect turbines have on the flow and may overestimate theresource most greatly in locations identified as being most suitable for genera-tion. Considering flow retardation as a result of the installation of turbines, gen-eration has a theoretical physical limit of 38% of the fluid power of the channel(0.38ρQ0∆Htide) corresponding in a reduction in bulk flow to 58% of natural.While the physical limits can be determined from channel and flow parameters,the extent of tidal development will likely be controlled by secondary constraintsrelated to environmental regulations or economics. The largest limit is related tothe permitted reduction in bulk flow. If only small reductions are permitted, thepotential for a channel to generate power is significantly reduced. This effect isamplified by any inefficiencies or extraneous drag associated with the turbines thatdiminish the resource without adding to generation. Furthermore flow reductionwill result in a five-fold decrease in per-generator output compared to the naturalconditions if considered at the physical limits of generation. This stepped flow re-duction could dramatically increase production costs and trigger disputes betweenstakeholders if development is incrementally staged.Tidal power clearly has the ability to generate some amount of electricity, how-ever the capacity for tidal power to produce at the utility-scale may be dramaticallyoverstated. Governments, developers, regulators and utilities must therefore becognizant of the physical, environmental and practical implications of this form of43generation when planning energy development. Policy and project planning shouldtherefore not be based on resource inventories utilizing inappropriate or oversim-plified methods such as the summation of kinetic energy flux.2.2 Impacts of Flow Reduction and Technical Learningon Levelized Cost of EnergyThe previous section discussed the physical and policy limits of extracting tidal en-ergy. The production of electricity at the large utility or industrial scale is, however,also contingent on the total cost of that energy being economically competitive.There are a number of factors that affect the cost of tidal energy (further dis-cussed below in Section 2.2.1) but is largely driven by the amount of energy aturbine produces, balanced against the capital cost to install that turbine. If a giventurbine produces less energy but the capital costs remain constant, the cost of en-ergy will go up. If the capital costs go down but the turbine produces the sameenergy, the cost of energy will go down.If there are large numbers of tidal energy turbines installed, both cost and out-put will change. The capital cost will drop as industry learns better ways of con-structing turbines. This is countered by drops in turbine product as the tidal flowsslow. The literature presents a model for technological learning (cost reduction),while a model of flow reduction was previously discussed in Section 2.1. Thesetwo models are combined to investigate their competing effects if large numbers oftidal turbines are constructed.442.2.1 Calculation of Levelized Cost of EnergyThe factors that affect the tidal power’s levelized cost of energy (LCOE) are thecapital cost of the turbines, the amount of energy produced, the cost of mainte-nance, the lifetime of the project and the interest rate payable on project financ-ing. LCOE is a useful parameter for assessing the economic feasibility of energyproduction methods as it represents the full price of energy once all the costs areconsidered. Using the notation of Short et al. (1995), LCOE is calculated as:LCOE =CapEx ·CRF +OpExAEP(2.19)Where CapEx is the overnight capital cost, OpEx is the annual operating ex-pense, AEP is the annual production of energy and CRF is the capital recovery fac-tor. The CRF represents the fraction of the project capital expense that is attributedto each year of operation, taking into account borrowing costs. It is defined as:CRF =i(1+ i)n[(1+ i)n]−1 (2.20)Where i is the cost of capital and n is the project lifetime. The CRF method isexactly equivalent to project analysis using discounted cash flows.Consider the recently proposed project by Marine Current Turbines, a UK tidalturbine manufacturer that is a considered a leading technology developer (WorldEnergy Council, 2013). The proposed array at Kyle Rhea is the installation offour 1.2MW SeaGen turbines at a cost of £35 million (at current exchange, ap-proximately $55M CAD). The company has stated that each of these turbines willproduce 6000MWh per year (Marine Current Turbines, 2013). Operating expense45for tidal stream turbines is estimated at £0.16M - £0.31M /MW/year ($250,000 -$485,000 /MW/year) (Ernst & Young LLP, 2010).As an illustrative example example, assume a 10 year equipment life and a 15%discount rate, the LCOE from one SeaGen turbine would be:CRF =i(1+ i)n[(1+ i)n]−1 =0.15(1+0.15)10[(1+0.15)10]−1=0.199LCOE =CapEx∗CRF +OpExAEP==13,750,000 ·0.199+300,0006000=$507MWhr2.2.2 Learning FactorsTh example LCOE of $507 / MWhr is approximately ten times the North Amer-ican spot price for electricity. The LCOE being above the market price is notunsurprising given tidal turbines are a new technology. The first units of many newtechnologies are often more expensive than the competition, but the price of theunits drops through volume and technological improvement.The relationship between number of installed capacity and unit price has beenwell studied with a number of empirical relationships being proposed (Henderson,1984; International Energy Agency, 2000; van den Wall Bake et al., 2009). Us-ing the methodology from the US Energy Information Administration (2013a), theunit cost of a technology will drop by some factor for each doubling of the in-stalled base. While EIA does not suggest figures for tidal power, the relationshipfor offshore wind is assumed to be the most relevant. It is important to note that46learning curves predict cost reductions during the development of a single technol-ogy. Marine Current Turbines’ device is a horizontal axis of conventional designand therefore these learning curves predict cost reductions for horizontal axis tur-bines. A shift in extraction technology to a different architecture could allow forcompletely different learning curves to be followed.For the first three doublings in capacity (up to 8 units), the price will reduce by20% for each doubling. For the next five doublings (up to 256 units) the price willreduce by 10% for each doubling. After 256 units are installed, the price will dropby only 1% for each doubling.Instead of a discrete cost reduction occurring at a doubling of capacity, costreduction is reformulated as a continuous function. As the machine enters massmanufacture, its unit cost drops, but levels-off at approximately 29% of its first-installation cost.If the calculation of LCOE from section 2.2.1 is redone with the CapEx andOpEx being reduced to 29% of original, the SeaGen turbine will produce electricityat $146/MWhr. While this rate is not competitive in the North American market, itis competitive in the European electricity markets.If some assumptions are made about an increase in cost of energy due toscarcity, tidal energy could be competitive in the future if the limitations of theresource are ignored.2.2.3 Interaction of Flow Reduction and LearningThere are two competing effects that are used in the model below to determine theLCOE of tidal power as a function of installed capacity:• The capital expense of the equipment, which drops as more turbines are built.47• The installation of multiple turbines reduces the velocity of the flow. Giventhat the amount of energy than can be produced is proportional to the cube ofthe flow velocity, flow reduction will reduce power production per turbine,leading to increased LCOE.MethodologyThe learning and cost reductions associated with total number of installed at theglobal scale, as technical learning is not confined by geography. Therefore, thesimulation flow retardation must also be global to capture the interplay of the twoeffects. There is, however, no global assessment of tidal energy resource. As hasbeen demonstrated, the inventory of kinetic energy flux is not relevant in assessingthe behaviour of the flow after turbines have been installed. As a result, to developa global resource, directed modelling of every suitable tidal channel would have tobe undertaken, therefore an assumption of the global resource must be made.Tarbotton and Larson (2006) identified Canada’s 50 largest tidal sites, whichwill be used as a exemplary set of the physical characteristics of developable sites.This set will be used as a basis of analysis of LCOE predictions. Tarbotton andLarson (2006) gives the maximum flood and ebb water velocity as well as the char-acteristic depth and breadth of the channel. What Tarbotton and Larson (2006)lacks however, is the pressure gradient (head difference) driving the flow throughthese channels. Without the head difference, it is not possible to calculate the avail-able resource of this set of channels. Given that the global resource is not known,the proceeding analysis will be conducted with the global tidal energy resource asan independent variable.The head difference will be calculated using the relationship for bottom friction48as given in Section 2.1 and recreated below ashI = fLRhV 22g (2.21)where V is the velocity in the channel and Rh is the hydraulic radius (Rh = A/P). Inorder to scale the sites, the friction coefficient, f , will be assumed to be 0.01 acrossall channels. The length of the channel, L, will be scaled geometrically with thechannel’s breadth. This scaling factor allows the global resource to be adjusted.The channels’ responses to adding turbines is simulated using a model writtenin Python. The simulation performs the following steps:1. One 20m turbine (the diameter of a SeaGen turbine) is added to the channelin the set that has the highest velocity.2. The resistance caused by the turbine is calculated and the velocity is reduced.Using the reduced velocity, the power generated is calculated and averagedover the a sinusoidal tidal cycle. The thrust and performance coefficients ofthe turbines are assumed to be 0.7 and 0.35 respectively.3. The new cumulative AEP from the turbines in all the channels is calculated.4. The incremental cost of the addition of this turbine is calculated using theprogress ratio method discussed above. This is added to the cumulative costof all installed turbines.5. The levelized cost of energy is calculated using the cumulative energy pro-duction and the cumulative installed cost.496. Repeat these steps until there are no channels in the set with a velocity above1.5m/s.These steps are repeated for varying magnitudes of global resource by varyingthe channel length.Results and DiscussionThe simulation of technological learning and flow reduction is computed for threemagnitudes of extracted tidal energy resource: 106, 107, and 108 MWhr/yr. Theselimits refer to the physical limits of η . Policy limits to flow reduction as discussedin the previous section are not imposed in this analysis. To reach maximum energyextraction in the 108 MWhr/yr scenario would require the installation of approxi-mately 13,000 tidal turbines. See figure 2.5.The extracted AEP differs from the global resource due to both the relation-ship between the thrust and performance coefficients, and the cut-in speed of theturbines. The coefficient of performance is assumed to be half of the coefficient ofthrust. This means that the turbine will impart resistance on the flow in addition tothe energy extracted from the flow. Physically, the difference represents frictionaldissipation and the drag on turbine support structures and ducting. Additionally, noturbine is installed in channels where the peak velocity is below 1.5m/s. Turbinesare not able to generate below a certain velocity, a cut-in speed. 1.5m/s is an oftencited value for this limit.The LCOE for the three conditions is shown in Figure 2.6 and 2.7. As ex-pected, the initial LCOE is well beyond market rates for all cases. The LCOEdrops as technical progress reduces the CapEx of the turbines. As more turbinesare installed, the velocity, and thus the kinetic energy and power output per tur-50108 MWhr/yr107 MWhr/yr106 MWhr/yrFigure 2.5: Annual Energy Production vs Installed Turbinesbine, begins to drop. The drop is less pronounced the larger the global resource.See Figure 2.8.The effects of the two competing forces yields a minimum LCOE for each sce-nario, which are shown in Table 2.2. As can be seen, if tidal energy must competeon a cost basis with market electricity rates, it is highly dependent on the size ofthe global resource.An important feature to note is the difference in capacity between the maxi-mum level of energy production that is physically possible, and the level at which51108 MWhr/yr107 MWhr/yr106 MWhr/yrLCOE (¢/kWhr)Figure 2.6: Levelized Cost of Energy vs Installed Turbines. The jaggednessof the curves is an artifact due as the model imposes an integer numberof turbines.minimum LCOE is found. The cheapest energy occurs when only a small factionof the resource has been extracted.Table 2.2: Minimum cost of energy for global resource levels.Global Resource Minimum LCOE106 MWhr/yr $218 / MWhr107 MWhr/yr $91 / MWhr108 MWhr/yr $62 / MWhr52108 MWhr/yr107 MWhr/yr106 MWhr/yrLCOE (¢/kWhr)Figure 2.7: Levelized Cost of Energy vs Annual Energy ProductionIf the tidal channels of the world can only support a few hundred turbines,then tidal energy is likely not going to become economically competitive. If manythousands of turbines can be constructed and installed, this form of energy maybecome economical. Until detailed modelling of the majority of tidal channels isundertaken however, the long-term feasibility of tidal energy technology will beuncertain.53108 MWhr/yr107 MWhr/yr106 MWhr/yrFigure 2.8: Annual Energy Production per Turbine vs Installed Turbines54Chapter 3Physical Modelling of TidalEnergy ResourceChapter 2.1 presented an analytical model to relate increased turbine resistanceto reduction in tidal flow. This reduction in tidal flow then lead to limitations ofpower generation. In order to validate the analytical model, a physical experimentwas constructed in the Hydraulics Laboratory of the UBC Department of CivilEngineering.Flow in a natural tidal channel is driven by a pressure gradient established bythe large tidal wave in the ocean. Turbines installed in this channel will developa head difference which takes away from the head driving the flow through thechannel, causing flow retardation. The power produced by the turbines is the prod-uct of the head developed and the flow rate. In this experiment, the tidal flow ina channel was simulated using actuated weirs to provide a head difference acrossa flume. Simulated turbine resistance is incrementally added using a shutter gate,55and the flow rate and surface elevations were recorded. These measured valueswere related to the analytical model parameters in order to validate the model.3.1 Boundary Conditions and the Conceptual ModelThe physics driving the flow of water in tidal channels is inherently different fromthat normally studied in the laboratory environment. As was discussed in the pre-vious chapter, the flow rate through a tidal channel is driven by a pressure gradientand restricted by friction and turbine resistance. This is contrasted with most lab-oratory flow channels or flumes that have flow produced by a pump which willbe near-constant. In order to simulate energy extraction in a tidal channel, a newexperimental setup must be developed. As such, the functional requirements of theexperimental setup were conceived as:• Geometry that reasonably represents a typical tidal channel found in nature;• Ability to have a continuously variable and controllable pressure gradient(head difference) across the length of the channel;• Provision to continuously vary the amount of resistance to flow in order tosimulate turbine energy extraction; and• Ability to measure the amount of energy extracted from the flow resistanceelement.A review of the literature did not reveal any similar experimental setup beingemployed at other institutions, therefore the conceptual design of the experimentwas from scratch.56It was determined early that using gravitational forcing would be difficult toreplicate in the laboratory setting due to scaling effects of gravitational attractionand lack of ceiling height. Conventional pumping is the practical alternative to de-livering water to the experiment. As such, the pressure gradient across the channelwould need to be directly controlled through either active or passive means.Inspiration was taken from a constant head tank, a long-used device to providea near-constant head to an experiment. See Figure 3.1. There are two caveats tothe use of constant head tanks: One is that the flow of water to the experiment mustbe significantly less than the inflow and drain from the constant head tank. Thesecond is that the head is not exactly constant. Changes in experiment flow thatwould otherwise change the head are counteracted by increased flow to the drain,which can only occur due to increased flow over the lip of the drain pipe (it actssimilarly to a weir). In many scenarios, the constant head tank is mounted suchthat this small change in tank level is insignificant compared to the total head, andcan be neglected. Additionally, the common pipe-based construction of a constanthead tank is fixed, such that only one tank level is achievable.A hybrid active/passive control system is utilized. At each end of the channelthere is both a water inflow and a water outflow. See Figure 3.2. In effect, each endof the channel behaves as a constant head tank. In steady state, the rate of waterinflow is exactly equal to the water outflow over the weir. Referring now to theweir at the right hand side of Figure 3.2, consider a demand for water flow towardsthe channel (to the left). The flow to the channel causes a deficit at the weir andthe water level will decreases slightly. Given that water flow over a weir increaseswith the water level above the weir (to the power of 1.5), an equilibrium will bereached between inflow, channel flow and drain flow. Provided the inflow and drain57TO EXPERIMENTTO DRAININFLOWFigure 3.1: Constant Head Tankflow are significantly larger than the channel flow, the change in water level will beslight. This forms the passive portion control system for water level.INFLOWINFLOWWEIRSIMULATEDTURBINESWEIRTHIS WEIR IS 90˚ TO FLOW.Figure 3.2: Conceptional model of water flume to be used.58Unfortunately due to the pumping flow rate requirements in the channel, it isnot practical to provide inflows and drain flows that are sufficiently larger than thechannel flows. As a result, an active control system must be developed to augmentthe passive control of water level. Accomplishing this requires actuation of thelevel of the weir. The water level at the end of the channel is measured and fedinto an electronic controller. The controller then changes the height of the weirsuch that a constant water level is maintained. These two control paradigms workin unison: the passive control reduced the higher frequency errors due to surfacewaves while the active controlled eliminates steady-state error.The effect of turbines on the flow must be imposed. A number of options wereconsidered including scale-models of turbines and porous discs. Scale models ofturbines were rejected due to the difficulty in controlling and varying the shafttorque, which affects the turbine’s resistance. It was also felt that sufficient resis-tance could not be generated in the space constraints of the chosen flume. Porousdiscs were rejected due to their lack of continuous variability and automatic con-trol.A shutter plate was selected as the ideal mechanism to provide flow resistanceas it is possible to vary the resistance from slight to nearly complete blockagethrough a continuously varying range. See Figure 3.3. In order to quantify theeffect of the resistance element, the head drop across the shutter plate must bemeasured as well as the flow rate through the channel.Flow rate measurement through the channel is by means of a two dimensionalintegration of a number of point measurements of velocity. Due to the large numberof data points required to develop the curves proposed in the theoretical model,and to improve accuracy and repeatability, these velocity measurements need to be59Figure 3.3: Rendering of the shutter plate design60automatically collected.Combining these sub-components with a controller that orchestrates the oper-ation gives an overall system that meets the requirements discussed above. Seefigure 3.4.INFLOW INFLOWWATER LEVELSENSORSHUTTERPLATEFLOWSTRAIGHTENERFLOWSTRAIGHTENERWATER LEVEL SENSORWEIRWEIRFLOW MEASUREMENT15mFigure 3.4: Arrangement of flume components. Bottom: scale of the flume.Top: detail of components.In order to conduct an experimental run, the operations below are performed.See figure 3.5.1. Initialize the experiment by starting all required equipment.2. Input the desired water levels at the boundaries.3. Run a control loop to adjust the downstream weir to attain the desired waterlevel.614. Run a control loop to adjust the upstream weir to attain the desired waterlevel.5. Repeat from 3 if the downstream water has changed.6. Measure head drop across the shutter plate.7. Measure the flow rate through the channel.8. Change the resistance provided by the shutter plate.9. Repeat from 3 until all shutter plate settings have been tested.3.2 Development of Physical Infrastructure3.2.1 Tank StructureTo implement the above control strategy, various components and sub-systemswere developed. The development process is discussed in the following sections.The flume to be used for the experiment is a 0.72m wide, 15m long water chan-nel installed in the Hydraulics Laboratory in the Rusty Hut building. See Figure3.6. In its original configuration, the flume is supplied with water from two cen-trifugal pumps, one 30HP and one 40HP, located in the basement of the building.The water is delivered to a head-pond at the upstream end of the channel via twovalved 10-inch (nominal) pipes and exit through a diffuser. The pumps draw wa-ter from a large underground sump, to which the discharge from the channel isreturned.The flume itself is built upon a pair of large steel channels that run the length.A hinge and a power screw allow for adjustment of the incline of the flume. For62InitializationSet Initial Desired Water LevelsAdjust Downstream WeirD/SWater Level Attained?Adjust Upstream WeirU/SWater Level Attained?D/SWater Level Changed?Measure Flow RateMeasure Head DropFinal Water Level?End NOYESNONONOYESYESNORecordDataNext Water LevelFigure 3.5: Process flow for experiment control system.63Figure 3.6: Two flumes in the Hydraulics Laboratory. The large yellow flumewas used for this experiment.this experiment, the flume was set to level. The flume structure, sump, upstreamwater outlet and pumps are the only components that were used in their existingconfiguration and were not built or modified for the purposes of the experiment.3.2.2 Water SupplyThe means of providing water to each end of the experiment was an effort in multi-variable optimization as available equipment, cost, and experimental flexibility areall factors. The design began with a performance model of the existing pumpinginfrastructure in Microsoft Excel, using the manufacturer’s pump curves and pub-lished values of loss coefficients. The initial values were validated against max-imum attainable flow in the unmodified flume and provided good general agree-64ment.The only feasible routing of the supply line for the downstream end of theflume was along the south wall of the Hydraulics Lab, then up and passing over theadjacent wave tank and walkway. The fittings and pipe lengths for this routing wereentered into the model, allowing calculation of water delivery rates for various pipesizes and pump configurations. Due to the availability of surplus pipe for the run, a6-inch nominal line was selected. The performance model indicates that with thisconfiguration, 0.127 m3/s can be delivered to the downstream end of the flume.The existing piping was cut and the new line connected in via a saddle flangeand secured to brackets bolted into the building structure. A wood structure sup-ports the PVC pipe as it passes over the wave tank and walkway before connectinginto a gate-valve and diffuser. See Figure 3.7c. The end of the pipe run is supportedby a 3/8-inch mild steel plate spanning between the gunnels of the flume structure.See Figure 3.8.Water injection into the upstream end of the tank was via the existing diffuserarrangement. See Figure 3.9.3.2.3 WeirsAn electronically actuated weir controls the outflow and water level at each endof the flume. The downstream weir was designed and constructed first and theupstream weir has a modified design, taking the learnings from the downstreamweir.65(a) Connection to LaboratorySupply(b) Pipe Run Along Wall(c) Pipe Spanning Wave Flume and Walkway, with SupportFigure 3.7: Water Supply to Apparatus66Figure 3.8: Downstream Valve and Steel Support PlateFigure 3.9: Upstream water inlet and diffuser located in upstream inflow tank67(a) Downstream weir from above (b) Downstream weir from outflowFigure 3.10: Downstream WeirDownstream WeirThe downstream weir is composed of three 3/16-inch acetal plates, two tracks,an actuator, and an aluminum frame. See Figure 3.10. The acetal plates slide inthe tracks (composed of UHMW plastic). The materials selected were to providea low coefficient of friction and robustness in a wet environment. An actuator,driven by an electric motor raises and lowers the weir. The design of the weirunderwent three iterations due to shortcomings of the initial design, primarily dueto difficulties with the actuator.The final iteration of the actuator was a commercial electric linear drive. Topower the linear actuator, an H-bridge and microcontroller were utilized. For posi-68tion feedback, a combination of a reflective optical sensor and linear potentiometerwere used. The robustness of this actuator came at the cost of performance as theactuator speed is quite slow. As a result, these actuators preclude the ability toconduct experiments with dynamically changing pressure gradients.Upstream WeirThe location of the upstream weir was slightly different than that of the downstreamweir due to the need to channel the spill water back into the sump. A hole was cutthough the 1/2-inch steel side wall of the head tank by department technicians. Itis through this hole that the discharge water flows and in it the weir is mounted.A wooden containment structure was built outside of the head tank to direct thedischarge water into a under-floor trough to return to the sump. See Figure 3.11.The size of the diversion structure was limited by the need to maintain walkwayaccess. Ideally the structure would have been larger as there was some choking andoverflow at very high discharge rates.The upstream weir mounts over the hole in the head tank. Instead of acetalsheet as the weir material, 3/8-inch marine plywood is used and reinforced withaluminum angle. The tracks in which the weir slides is made from aluminiumchannel. This configuration provides improved stiffness over the acetal weir andhas fewer issues of jamming. The upstream weir is driven by the same electriclinear actuator as was used for the downstream weir. See Figure 3.12.3.2.4 Flow Resistance ElementThe flow resistance element is the component that simulates the energy extractionfrom tidal turbines. It must be continuously variable and provide a wide range69Figure 3.11: Closure built for upstream outflow. Leads to in-floor trough tobe returned to the sump.70(a) Looking into head tank from abovecontainment structure(b) Looking into containment structurefrom above head tankFigure 3.12: Upstream weirof resistances from nearly uninhibited to completely blocked. A shutter plate wasselected as the resistance element. See Figure 3.13.The shutter plate consists of six slats constructed by machining 1x4 lumberstock. These slats must be thin in order to provide little resistance when fullyopen but strong enough to support the hydrostatic pressure across them when fullyclosed. The slats are painted to increase their lifespan while in contact with water.A 1/4-inch stainless steel rod is inserted and epoxied into each end of each slat.This rod acts as an axis about which the slats rotate. The rods are mounted intoPTFE bushings secured to the bottom of the flume and to a beam spanning betweenthe gunnels of the flume.71Figure 3.13: Shutter plate slats.At the top mount point, the rods extend through the bushings and are exposedfor control. A cam arm is affixed to each rod to pivot the slat. All the cam armsare connected through linkage bars so that all the slats rotate in unison. See Figure3.14. A 1/4-20UNC threaded mount is affixed to the linkage arms and acceptsa lead screw. This lead screw is driven by a stepper motor (reclaimed from thefirst iteration of the weir) allowing automatic control of the magnitude of flowresistance. Limit switches ensure that the absolute position of the slats can bedetermined by the motor drives.The motor communicates with the control computer via RS-485 serial whichis converted to USB via an adapter.3.2.5 2D ProfilerTo determine the amount of energy being dissipated by the shutter plate, the flowrate through the channel must be determined. This is via 2D integration of veloc-72Figure 3.14: Linkages that actuate the slats of the shutter plate.ity measurements over the cross section of the channel. The measurement of thevelocities should not meaningfully change the resistance of the channel.Existing in the inventory of the Hydraulics Laboratory is a two-axis traverse.This device combines linear slides with ACME lead screws and stepper motors tomove a plate in X and Y directions. See Figure 3.15. The drive electronics for thisdevice were built to accept pulse inputs from a parallel port - a peripheral connec-tion not found on modern computers. As a result, a new electronic interface neededto be developed. These stepper motor drives accept a 5V pulse and direction input;each pulse input causes the motor to turn 1/200 of a rotation. As a result, a highspeed pulse train is needed to attain any meaningful speed. An Arduino microcon-troller was used to interface the control PC (via USB) and provide the pulse anddirection inputs to the motor drives. Code was written for the microcontroller thatdealt with motor acceleration, deceleration, positioning, homing and limits.To measure the water velocity, a Sontec acoustic doppler velocimeter (ADV)73was affixed to the traverse. This device provides X, Y, and Z velocities in a 1cm3measurement volume. The ADV reports its data back to the control PC via anRS-232 serial interface.Figure 3.15: 2D Traverse that supports the ADV.74To calculate the flow rate, the control PC commands the traverse to a knownposition in the channel cross section. Velocity measurements are taken at this pointand stored with their spatial coordinates. The traverse is then commanded to thenext position, and so on, until an array of velocities is collected. This array canthen be spatially integrated. See Section 3.3.3.2.6 Water Level Elevation MeasurementIn both the calculation of energy dissipation and controlling the pressure gradientat the boundaries, the elevation of the water level must be measured. A numberof devices were considered and an ultrasonic distance sensor was selected. Thesedevices emit an ultrasonic chirp, which, when reflected off objects or interfacesis later received by the sensor. By timing the interval between the emission of thechirp and the reception of the reflection, the distance to the object can be calculated.These sensors are pointed at the water surface, which is the source of the soundreflection. This allows a non-contact measurement of the water surface relative toabsolute elevation.Four SensComp model Mini-A ultrasonic sensors are utilized in the experi-ment: one at each of the upstream and downstream ends of the flume, and one oneach side of, and immediately adjacent to, the shutter plate. The sensor is mountedon a rod, held by a support that spans the gunnels of the channel. A custom housingand power supply circuit was built for each of the sensors. See Figure 3.16.The ultrasonic sensors output a 5V analog signal, proportional to the distancemeasured. This signal is carried via shielded cable to a data acquisition system atthe control PC.75Figure 3.16: Ultrasonic sensor near downstream flow conditioning grid.763.2.7 Bottom ResistanceThe bottom and sides of the channel are comprised of PVC and acrylic sheets withsmooth surfaces. Left unaltered, a head difference across the channel of only 1-2mm would consume all of the available flow from the inlet pipes. Such a conditionwould produce unacceptable measurement errors as the resolution of the ultrasonicsensors is as large as 0.3mm. Increasing the bottom friction means a large pressuregradient results in a smaller (and attainable) flow rate.Bottom friction in the flume is enhanced by the installation of a lattice structureto produce a large effective roughness coefficient. The lattice is constructed byfastening together pieces of 1.5-inch by 0.75-inch dimensional lumber (1x2s) onedge. See Figure 3.17. The transverse members are spaced at 0.3-m intervals.The lattice runs the length of the channel and is secured to the bottom to preventflotation.The enhanced bottom friction allows for moderate differences in head to beestablished without exhausting the available inlet flow.3.2.8 ElectronicsDue to the large number of subsystems, power and data runs require organization.The electrical system is centred on a control bench utilized to run the experiment.The bench is fed with a 120VAC power line from the building mains. The 120VACis received by a panel installed in the bench. Power is available via sockets on thebench for computers and other hardware. Power being sent to electrical systemsaway from the bench is hardwired through accessible switches. As a result of thislayout, power can be shut off to any component by the operator at the bench. Thisprovided both safety and operational advantages, as malfunctioning components77Figure 3.17: Friction lattice on bottom of channel. Looking upstream fromdownstream weir.78can be easily de-powered (preventing damage) without the need to shut down othersystems like the control PC. See Figure 3.18.DATA COLLECTION PCBenchBox2D TRAVERSE DRIVESWITCH 1SWITCH 2SWITCH 3DOWNSTREAM WEIR P/SULTRASONIC SENSOR P/SWEIR ACTUATORSHUTTER PLATE P/SULTRASONIC SENSORS P/SSTEPPER MOTORSONIC SENSORSONIC SENSORSONIC SENSORDOWNSTREAM WEIR P/SULTRASONIC SENSOR P/SWEIR ACTUATORSONIC SENSORMAIN POWER SWITCHADV CONTROLLERLAB MAINSNETWORK ADAPTER P/S12VDC18VDC18VDC40VDC18VDC12VDCAll lines120VAC unless noted.Figure 3.18: Power connections for the apparatus.The Hydraulics Lab is an electromagnetically noisy environment, so digitalcommunication is the preferred transmission mode. To the extent possible, Univer-sal Serial Bus (USB) has been used due to its scalability and robustness. Due tothe relatively long transmission distances between the control PC and the micro-79controllers at the weirs, repeaters are required. See Figure 3.19.Analog data from the ultrasonic sensors is transmitted over shielded 2 conduc-tor wire. BNC connectors allowed for easy mating between the sensors and theMeasurement Computing PCI data acquisition card.DATA COLLECTION PCDAQ CARDUSB HUBPCIUSBSONIC SENS 1SONIC SENS 2SONIC SENS 3SONIC SENS 42D TRAVERSE MICROCONTROLLERMDRIVE STEPPERU/S MICROCONTROLLERADV CONTROLLER0-5V AN0-5V AN0-5V AN0-5V ANRS-232RS-485 CONVERTERD/S MICROCONTROLLERUSBUSBRS-485USBUSBH-BRIDGEH-BRIDGESTEPPER DRIVE5V DIG5V DIGENCODERENCODER5V DIG5V DIG5V DIGFigure 3.19: Data connections for the apparatus.The control software and process is discussed in Section 3.3.803.3 MethodsThe purpose of this experiment is to test the validity of the analytical model pre-sented in 2. Thus, the specific physical data collected represents the variables ofthe model. This are notably:• ∆Htide - head difference between the ends of the channel• Q - the flow rate in the channel• Q0 - the flow rate before any turbine resistance is applied• hI - head lost to bottom friction• hT - head drop across the turbines• kI - bottom friction coefficient• kT - turbine resistance coefficientUsing these parameters η , the extraction efficiency, can be calculated such thatthe analytical and physical models can be compared. The primary independentvariable is the turbine resistance coefficient, kT . kT is varied through actuation ofthe shutter plate. One experimental run is comprised of a set of kT while holdingother independent variables constant.The other independent variables are ∆Htide and kI . ∆Htide can be set in thecontrol software (subject to the limits of the pumping capacity). kI is fixed in thissetup, but may be varied by changing the lattice structure on the bottom of theflume. hT , hI , Q and Q0 are dependent variables.81To calculate these parameters, many measurements are taken, combined andprocessed. Table 3.1 shows what physical measurements are used to calculate theanalytical parameters.Table 3.1: Experimental measurements used to compute analytical modelvariables.Model Parameter Physical Measurements∆Htide Ultrasonic measurements at each end of the flumeand initial arrangement measurements.Q ADV velocity and 2D profiler positionQ0 ADV velocity and 2D profiler positionhI Ultrasonic measurements at each end of the flumeand initial arrangement measurements.hT Ultrasonic measurements on each side of the shut-ter plate, initial arrangement measurements, andQ.kI hI and Q, as described in Chapter 2.kT hT and Q, as described in Chapter 2.3.3.1 Calculation and Control of ∆Htide∆Htide is the different in elevation of the water surface between the water adjacentto the upstream and downstream weirs. Its is important that the water surfacemeasurement be taken sufficiently far from the weir such that the measurement isnot influenced by the hydraulic control at the weir crest. The measurement takenduring the experiment is the distance from the ultrasonic sensor to the water surfacewhereas the measurement needed is the elevation of the water surface relative to anabsolute and universal datum. This conversion is performed using measurementsand calibrations taken before the experiment was started as follows.1. Using a transit, measure the elevation of the floor of the flume to the datum,82Zdatum→ f loor.2. Using two boxes, each of a known height and with an acoustically reflectivesurface, perform a two point calibration of the ultrasonic sensor reading tothe distance from the flume floor. In this arrangement, the signal from thesensor will be reversed (negative calibration coefficient).From the calibration step, the raw sensor output can be converted from 0-5VDCto distance from the floor to surface, Z f loor→sur f ace. The water surface elevation isthen given byZsur f ace = Zdatum→ f loor +Z f loor→sur f ace (3.1)With the elevation relative to the datum, ∆Htide can be calculated as:∆Htide = Zupstreamsur f ace −Zdownstreamsur f ace (3.2)The sensor signal itself requires significant processing before is can be con-sidered useful. If the sensor emits a chirp and it does not receive a strong enoughreflection, its voltage output goes to 5V leaving artificial spikes in the data. Analgorithm was written to recognize these spikes and remove them.Due to the long distance of the cable run between the sensors and data acquisi-tion card in the control PC, there was significant electrical noise on the signal line.This is removed using a low-pass filter utilizing a Fast Fourier Transform (FFT).Surface waves and ripples are also considered noise and are averaged out througha 10-second measurement at 10 samples/second.When the experiment is started, ∆Htide is an input parameter to the controlsystem as both Zupstreamsur f ace and Zdownstreamsur f ace . At each control loop, the system actuates83the weirs to maintain these water surface elevations within the 0.5mm of their setpoints. This is accomplished using a P-controller algorithm.3.3.2 Calculation of QThe flow rate Q is calculated from multiple measurements of velocity collected bythe ADV and integrated as described in Section 3.2.5.Matlab code was written to parse the ADV serial stream, which includes XYZvelocities as well as XYZ correlation coefficients. The correlation coefficients rep-resent the strength of the acoustic return for that measurement. Measurements withcoefficients below 80% were dropped. The signal was then subjected to a FFT low-pass filter.The set of velocities at their corresponding XY coordinates is converted into abulk flow rate. To accomplish this, the cross section of the flume is divided intorectangles, each of which is applied to each velocity measurement. See figure 3.20.The flow rate through a single rectangle is the product of the associated measuredvelocity and the rectangle’s area. Bulk flow rate is the sum of the individual flowrates associated with each rectangle.The borders of rectangles are the bottom or sides of the flume, the free surface,or equidistant between the measurement locations.3.3.3 Calculation of hT and kThT is calculated using the two ultrasonic surface elevation measurements madeimmediately adjacent to the shutter plate, ZSP1 and ZSP2. These measurements arecalibrated and converted in the same manner as the weir sensors as described inSection 3.3.1.84VelocityMeasurementPointArea applied tomeasurementFigure 3.20: Cross sectional area applied to each velocity measurement forthe purposes of bulk flow calculation.hT is not simply the difference between ZSP1 and ZSP2 because the shutters de-velop a head difference even whilst fully open in a manner akin to natural channelresistance. This initial head difference should be subtracted from hT in the cal-culation of power, however it also varies as a function of flow. Therefore, hT iscalculated ashT = hmeasuredT − kSP0 ·Q2 (3.3)where hmeasuredT is the post-calibration head difference directly observed by theultrasonic sensors (ZSP1−ZSP2). kSP0 is the open resistance coefficient of the shut-ter, calculated as:kSP0 =hT 0Q20(3.4)where hT 0 is the measured head difference when the shutter is fully open and Q0 isthe flow rate during the measurement of hT 0.85kT is calculated askT =hTQ2=(ZSP1−ZSP2)− kSP0 ·Q2Q2(3.5)This presentation of kT relates head with flow, instead of velocity as presented inChapter 2, a variation of a factor of A2. This has no effect on the result however,as the independent variable being tested is kT/ki. ki is similarly formulated withrespect to flow and thus the factor of A2 cancels out. The choice to present withrespect to flow is to place emphasis on the measurement of power production asthe product of flow and head, given below.3.3.4 Calculation of hI and kIhI is calculated from the residual of ∆Htide and hT using the measurement tech-niques discussed above.hI = ∆Htide−hT (3.6)kI again the relationship between hi and Q askI =hIQ2(3.7)Note as above, this formation of kI is with respect of flow, not velocity.3.3.5 Power and ηThe power dissipated by the the shutter plate, the turbine analog, is fluid power lostas given by measured variablesP = ρghT Q (3.8)86This now allows calculation of η asη =Pρg∆HtideQ0=hT Q∆HtideQ0(3.9)3.3.6 OperationAll data are collected and logged automatically in the central control system, whichconsists of numerous Matlab scripts, functions and toolboxes running on a Win-dows PC with data acquisition hardware.An experimental run is conducted as follows:1. Calibrate the ultrasonic sensors.2. Set the shutter to full open (no generation)3. Measure the boundary water levels and actuate the weirs and iterate to attainthe water level set-points. Iterate until the set-points are stable.4. Measure the head difference developed across the shutter.5. Measure the flow rate through the channel.6. Incrementally close the shutter (simulating the installation of additional tur-bines).7. Repeat from Step 3.All the data collected is stored as Matlab binary files for post processing.873.4 Results and DiscussionOf the two independent variables, a matrix was formed from 3 ∆Htide values, eachwith between 10 and 16 shutter angles. One setting of ∆Htide is considered anexperimental run. The values of ∆Htide are shown in Table 3.2. Data from all theexperimental runs is shown at the end of this chapter in Table 3.3.Table 3.2: Surface elevations and ∆Htide for conducted experimental runs.Experimental Run Zupstreamsur f ace Zdownstreamsur f ace ∆Htide1 0.5306m 0.5065m 0.0241m2 0.5335m 0.5099m 0.0236m3 0.5460m 0.5230m 0.0230mFor each of these runs, the shutter plate was incrementally closed. This causedflow rate to decrease and the water levels to change. The response in water sur-face elevation for experimental run 3, a typical run, is shown in Figure 3.21. Theexperimental results are presented in numerical form in Appendix ??.0 10 20 30 40 50 600.520.530.540.55Shutter Angle (degrees)Water Surface Elevation (m)Upstream of shutterDownstream of shutterUpstream boundaryDownstream boundary ∆H tideABCDEFFigure 3.21: Change in water surface elevation along the flume as a result ofshutter closure during experimental run 3.Consider the water levels at a zero shutter angle, denoted in Figure 3.21 by88points A (adjacent to upstream weir), B (adjacent to shutter plate on the upstreamside), C (adjacent to shutter plate on the downstream side) and D (adjacent to down-stream weir). The head difference A→ B is due to bottom friction in the upstreamsection of the flume. The drop from B→C is the difference over the shutter platewhile open. B→C represents hT,measured as discussed in Section 3.3.3. The headdifference from C→ D is equivalent to A→ B but at the downstream end of theflume.As the shutter plate is closed, the water level differences between the weirs andthe near-side of the shutter plate decreases to the point of maximum closure, whereis no measurable bottom friction head loss. This can be seen at points E and F inFigure 3.21. Conversely, the head drop across the shutter plate steadily increasesas the shutter is closed. The water level between points A→ E and D→ F areenforced by the control system via the weirs such that A→D and E→ F are equalto ∆Htide.hT is calculated as the difference between the water surfaces immediately up-stream and downstream of the shutters and corrected. hT rises, asymptotically ap-proaching ∆Htide as the shutter produces full blockage. Concurrently, the closureof the shutter reduces the flow velocity, which reduces the bottom friction headloss as discussed above. See Figure 3.22a. These two behaviours are intuitivelyexpected in this flow regime.Dissipated power is the product of hT , Q, gravity and density, and is calculatedover the range of shutter positions. See Figure 3.22b. For low shutter angles, poweris limited by a small hT while Q has not been substantially reduced. Power reachesa maximum and falls for high shutter angles as Q is diminished.η , which represents the ratio of extracted tidal power to natural fluid power of890 10 20 30 40 50 6000.0050.010.0150.020.0250.03Turbine Head, h T(m)5075100125Flow Rate (L/s)(a) Shutter Angle (degrees)QhT0 10 20 30 40 50 60051015(b) Shutter Angle (degrees)Extracted Power (W)Figure 3.22: (a). Change in turbine head and flow rate as a result of shutterclosure. (b). Extracted power from flow as function of shutter angle.Data are from a Experimental Run 3.the system is is plotted as a function of Q in Figure 3.23. These data show goodagreement with the analytical model (Equation 2.7) and shown as a dashed line.The conducted experiments found maximal values of η in the range of 0.41-0.45,depending on the experimental run. This shows good agreement with the analyt-ical maximal η of 0.38. The 6-15% difference between analytical and measuredvalues for η are within experimental uncertainty, which will be further discussedin Section 3.4.1.Another important aspect in generating large amounts of electricity from tides900.30.40.50.60.70.80.9100.10.20.30.40.5Extraction Ratio ηFlow rate as a ratio of natural (Q/Q0)Predicted relationship fromanalytical model. Run 1Run 2Run 395% of naturalFigure 3.23: Extraction efficiency (η) versus flow rate. Each set of symbolsrepresents an experimental run.is the necessary changes in flow needed to realize this maximum power output. Theanalytical model predicts that the bulk flow rate will be reduced to 58% of natural,at the point of maximum production (see equation 2.10). Again, the experimentalresults show agreement. The maximal extraction peak occurs between 61% and64% of natural flow, depending on the experimental run.Reducing the flow in a marine system to roughly 60% of natural will likely bea unpalatable proposition in many jurisdictions. Therefore assessing what is therealizable resource is assessing the effects of limiting flow reduction. Hagermanet al. (2005) suggests that the channel’s flow must be within 95% of natural onaccount of environmental issues. If this condition is imposed on the experimentalresults, the maximum interpolated η is between 0.15 and 0.20. These values arelarger, but within experimental uncertainty of the analytical prediction of η = 0.09.913.4.1 Experimental UncertaintyAn uncertainty analysis will be conducted with a Monte Carlo simulation using abeta distribution of possible values for each measurement. This method is basedon Project Evaluation and Review Technique (PERT) analysis. Each measurementis assigned a best-guess, minimum value and maximum value as well as a shapefactor, α (set as 2.0).The Beta distribution is a two-parameter, bounded statistical distribution withthe mode of the PDF set as the best-guess. A Beta PDF is shown in Figure 3.24.3 3.5 4 4.5 5 5.5 600.10.20.30.40.50.60.70.80.91Figure 3.24: Beta distribution for a best guess of 4, minimum value of 3 anda maximum value of 6.92Surface ElevationSurface elevation measurements are made using ultrasonic sensors, which operateby reflecting pulses of sound off the water-air boundary and timing the return-time;this gives the relative distance between the sensor and the surface. The sensorshave stated accuracy of 1% and a resolution of approximately 0.3 mm but canhave a zero-drift of up to 1.5 mm as a function of air temperature and humidity;the position of the sensor is therefore calibrated off the floor of the flume beforeeach trial. The calibration is done with an empty flume, changes in temperatureand humidity due to the moving water will effect the speed of sound. The largestuncertainty is the relative elevation of the flume floor between sensors, as measuredusing a surveying level. Due to the arrangement of the laboratory and the opticsof both the flume walls and the level, it was not possible to resolve inter-sensorelevation beyond 1mm. This means that while the initial value of ht or ∆Htidemay be relatively inaccurate, we may be relatively confident in the magnitude ofchanges in ht or ∆Htide. For the purpose of analysis, the best-guess is the measuredvalue, minimum is 99% of measured (sensor accuracy) less 1.5 mm (0.3 mm forresolution, 1.2 mm for relative elevation and zero-drift), and maximum is 101% ofmeasured plus 1.5 mm. These ranges are used with the understanding that relativechanges will have better agreement.Flow RateVolumetric flow rate is a difficult parameter to accurately measure due to the non-uniformity of the velocity profile across the flume. For each flow rate measurement,25-36 velocity measurements are collected and averaged. The variability of veloc-ity profiles is shown in Figure 3.25. Temporally separated point measurements93required significant interpolation but there are trade-offs in experimental design.Increasing the number of points means there is significant time between the firstand last measurements. Reducing the length of each measurement would shortenthe total time, but induces the risk of bias from turbulence or surface waves. Theseruns use moderate measurement times (10-20 seconds each) with relatively widespacing, while still requiring 8-10 hours of run-time. For each measurement, thebest guess is the integrated velocity using all measurements, the minimum is 90%of the minimum specific flow vertical profile and maximum is 110% of the maxi-mum profile.Aggregate Data with UncertaintyThe data is parameterized using η as above and kT/kI being the resistance coeffi-cients for the shutter gate and channel, respectively.The aggregate results with the upper and lower quartiles are shown in Figure3.26. The analytical predictions fall well within the envelope of the experimentalresults. Through inspection of the results of the uncertainty analysis, the followingis concluded:• The relatively tight grouping of data is the result of the stability and relativeprecision of the sensors.• Virtually all data lies above the analytical prediction, yet follows the formof the analytical result very well. This indicates that the dominant source ofexperimental could be a static offset in one or more inputs.• As a percentage of measurement, there is greater uncertainty with flow ratethan water surface elevation. As a result, confidence bounds are much wider940 0.2 0.4 0.6 0.8 100.050.10.150.20.250.30.350.40.450.5Normalized Velocity (u/umax)Depth (m)Figure 3.25: Variation in the shape of the velocity profile over experimentsand positions.for small values of kT/kI where Q plays a relatively more significant rolethan hT in determining power.3.5 Conclusions and RecommendationsA new physical modelling apparatus was constructed to validate analytical mod-els for tidal power generation. The development of this infrastructure appears tobe unique and simulates boundary conditions that do not seem to have been mod-elled before. The engineering effort to construct and commission the setup was950 2 4 6 8 1000.10.20.30.40.50.6Extraction Ratio ηRelative Resistance kT/kIResults predicted by the analytical model.Figure 3.26: Fraction of total fluid power extracted (η) versus turbine resis-tance relative to natural channel resistance. Error bars represent upperand lower quartile of Monte Carlo uncertainty results. Each set ofsymbols represents an experimental run.significant and several subcomponents required multiple iterations.Significant emphasis was placed on developing the apparatus for the lowestcost possible. There were many examples where this was false-economy as thefrugality significantly delayed the overall research program.Several experimental runs were conducted and the results of the analyticalmodel are consistent with the physical model results within experimental uncer-tainty.96The analytical model indicates that the limit of 38% of the unimpeded fluidpower of the system may be extractable. The physical experiments produced amaximum of 41-45% of the natural fluid power, representing a 6-15% deviationbetween the models.Experimental data show that maximum production occurs when flow is reducedto 61-64% of natural. If bulk flow rates are held by regulation to near natural rates,this regulation will act as the control on the tidal power resource. The smallerthe reduction that is permitted, the larger the fraction of the physical resource thatcannot be realized. If flow reductions are kept small, extractible power falls quicklywith Q/Q0. These limits are highly sensitive to small changes in resistance and itmay be difficult to accurately predict the flow reduction until full-scale devices areinstalled.The experimental results show relatively good agreement with the analyticalpredictions for flow reduction.Additional work may be warranted in testing additional values for ∆Htide andvarious water depths, as this independent variable is not extensively explored.Furthermore, recalibration and duplication of the existing runs could illuminatethe uniform discrepancy between the model and physical measurements.97Table 3.3: Selected data from laboratory experiments.Run Zupstreamsur f ace ZSP1 ZSP2 Zdownstreamsur f ace Flow hT ki kT Power ηm m m m m3/s m - - W -1 0.531 0.5‘23 0.513 0.506 0.0754 0 4.397 0.000 0 01 0.531 0.523 0.513 0.506 0.0754 0 4.397 0.000 0.241 0.0141 0.531 0.524 0.512 0.506 0.0808 0 4.397 0.000 -0.099 -0.0051 0.531 0.525 0.512 0.506 0.0765 0.002 4.397 0.342 1.793 0.0981 0.531 0.526 0.511 0.507 0.0708 0.006 4.222 1.197 4.12 0.2291 0.531 0.527 0.509 0.506 0.0687 0.009 4.397 1.907 5.906 0.3121 0.531 0.528 0.51 0.507 0.0619 0.011 4.222 2.871 6.594 0.3691 0.531 0.528 0.509 0.506 0.0566 0.013 4.397 4.058 7.056 0.3971 0.531 0.529 0.509 0.507 0.0505 0.016 4.222 6.274 7.667 0.4331 0.531 0.53 0.508 0.507 0.0465 0.017 4.222 7.862 7.855 0.4351 0.531 0.53 0.508 0.506 0.0393 0.02 4.397 12.949 7.521 0.4131 0.531 0.53 0.507 0.506 0.0379 0.02 4.397 13.924 7.451 0.411 0.531 0.531 0.508 0.507 0.0321 0.021 4.222 20.380 6.554 0.3692 0.533 0.526 0.515 0.509 0.0772 0 4.027 0.000 0 02 0.534 0.526 0.515 0.509 0.0781 0 4.195 0.000 -0.155 -0.0082 0.534 0.527 0.515 0.509 0.0781 0.001 4.195 0.164 0.527 0.0282 0.534 0.527 0.514 0.509 0.0759 0.003 4.195 0.521 2.269 0.1212 0.534 0.529 0.513 0.509 0.0686 0.007 4.195 1.487 4.717 0.2572 0.533 0.529 0.512 0.509 0.0646 0.01 4.027 2.396 6.174 0.3352 0.533 0.53 0.511 0.509 0.0593 0.012 4.027 3.412 6.982 0.382 0.533 0.531 0.51 0.509 0.0534 0.015 4.027 5.260 7.857 0.4232 0.533 0.531 0.509 0.508 0.0482 0.017 4.195 7.317 8.205 0.4322 0.533 0.532 0.509 0.509 0.047 0.018 4.027 8.148 8.432 0.4482 0.533 0.532 0.509 0.509 0.0417 0.02 4.027 11.502 7.994 0.4292 0.534 0.533 0.509 0.509 0.0346 0.021 4.195 17.542 7.114 0.3842 0.534 0.533 0.509 0.509 0.0326 0.022 4.195 20.701 6.987 0.3712 0.533 0.533 0.509 0.509 0.0275 0.022 4.027 29.091 6.002 0.3252 0.534 0.533 0.509 0.509 0.0247 0.023 4.195 37.699 5.555 0.2943 0.547 0.538 0.528 0.523 0.1315 0 1.388 0.000 0 03 0.547 0.538 0.528 0.523 0.1307 0.001 1.388 0.059 0.654 0.0213 0.547 0.54 0.527 0.522 0.122 0.005 1.446 0.336 6.366 0.1983 0.547 0.542 0.526 0.523 0.0996 0.01 1.388 1.008 10.204 0.333 0.547 0.543 0.524 0.522 0.0849 0.016 1.446 2.220 13.067 0.4093 0.547 0.545 0.523 0.522 0.0681 0.019 1.446 4.097 12.818 0.4023 0.547 0.546 0.523 0.523 0.054 0.021 1.388 7.202 11.154 0.3583 0.547 0.546 0.523 0.523 0.0419 0.023 1.388 13.101 9.325 0.398Chapter 4Numerical Modelling ofHydrokinetic Generation Impactsnear Traditional Hydropower1As was discussed in Chapter 1, hydrokinetic generation refers to the use of free-stream turbines for generating power from directly from the moving water of rivers,without the construction of dams (Federal Energy Regulatory Comission, 2015).These turbines operate in the same manner as tidal turbines in that they convert thekinetic energy of the water to electricity, but the installation in rivers means thephysical drivers of the flow are markedly different from tidal channels.This chapter will discuss some implications of hydrokinetic power develop-ment in a river near an existing hydropower facility.1Note: The text of this chapter is adapted from a technical report submitted to and accepted bythe BC Hydro and Power Authority. It was specifically a study of a proposed installation of turbinesnear a conventional hydroelectric generation facility for the purpose of determining if BC Hydro’sshould allow the project. Some information has been excluded to protect confidential information.99Hydrokinetic power is electricity generated by extracting the kinetic energy ofthe water in a river as it moves through a turbine placed in the free stream. This iscontrasted to traditional hydropower, in that the necessary pressure gradient acrossthe turbine is developed by the resistance of the turbine and not imposed througha dam and penstock, or other structure. Presently, hydrokinetic turbines are pre-commercial, and a number of developers hope to test their devices in order to proveperformance (US Energy Information Administration, 2012)An attractive location to test a hydrokinetic device is near the infrastructure ofa traditional hydropower dam. These locations typically offer relatively constantand debris-free flows with easy connection to the electricity grid. Depending on thescale and location of the turbines, they may act either to augment the hydropowergeneration, or to detract from it. As a result, the impacts of these turbines shouldbe understood before a project is undertaken.As a case study, a hydrokinetic technology developer proposed such an installa-tion near BC Hydro’s Seton Generation Station. The proposed installation locationis Seton Canal, a channel that conveys water from Seton Reservoir to the intake ofSeton Generation Station. BC Hydro has concern that this project would negativelyimpact the conventional hydro-electricity production at this powerhouse. The re-mainder of this chapter will evaluate the impacts of hydrokinetic generation on atraditional hydropower facility.1004.1 Background4.1.1 Existing Seton Reservoir GenerationThe Seton Generation Station (SGS) is located approximately one kilometre southof Lillooet, British Columbia (50.67◦N, 121.92◦W) and is part of the Bridge RiverSystem. This station has an average output of 48 MW, generated from a 100 m3/sflow from Seton Reservoir, discharging into the Fraser River (Percel, 2010).Water is conveyed to SGS from Seton Reservoir by means of a 3.5 km canal,originating at the Seton Dam. A portion of the reservoir outflow is channelled intothe canal, while the remaining outflow is spilled into the Seton River which flowsinto the Fraser River. See Figure4.1. The proportion of flow into Seton Canal isnot controlled directly, but is a function of the amount of spill and reservoir level.SGSSetonReservoir Seton RiverSeton CanalLillooetFraser RiverFigure 4.1: Seton Reservoir, Canal and Generation Station. Image adaptedfrom Google EarthApproximately 0.3m of head is lost in conveying water from Seton Reservoirto SGS, corresponding to a frictional dissipation of 300-400 kW. Mean velocities101in Seton canal range between 1.1 and 1.9 m/s in a trapezoidal channel roughly 4-6m deep and 30 m wide.At steady state, the total energy flux into any system must equal the flux out;the form of these energy fluxes may change, however. The Seton Canal systemhas three primary out-fluxes, electricity from SGS production, kinetic energy inthe discharged water (a function of SGS output and configuration) and frictionaldissipation in the canal and SGS. Given a constant energy in-flux, the electricityextracted by hydrokinetic turbines must therefore be a substitution of one of theforms of energy out-flux.There would be no net production benefit of the project if the electricity produc-tion from hydrokinetic turbines was less than the amount their installation reducesSGS production. Thus, the increase in generation should be derived largely fromreduced frictional dissipation if positive net-power is to be achieved.4.1.2 Theoretical BasisThe effect of a disturbance on the flow is a function of the flow criticality. If re-sistance (such as that from a turbine) is installed in subcritical flow, the resistancewill cause back-watering upstream (downstream controlled). If the resistance isinstalled in supercritical flow, the impact will only be felt downstream (upstreamcontrolled). The flow in Seton Canal is subcritical and therefore downstream con-trolled. The head difference developed by the hydrokinetic turbine will propagateupstream until one of the following conditions is encountered:a. a hydraulic jump, above which the flow is supercritical and the depth and ve-locity is independent of downstream conditions. See Figure 4.2a. This typeof flow occurs in steep channels or downstream of narrow constrictions or102weirs. On site visits and aerial photography of Seton Canal indicated that theentire canal appeared to be subcritical flow.b. the flow reaches uniform depth, a condition that occurs when the gravity driv-ing the flow is in equilibrium with friction that retards the flow, resulting inconstant depth. See Figure 4.2bConsidering only changes in depth, if the backwatering does not encountersupercritical or uniform flow and instead encounters an externally controlled waterlevel, a discontinuity would occur (See Figure 4.2c). A controlled water level existsat Seton Lake, therefore both the depth and the bulk flow rate in the canal wouldchange if the flow is not controlled by supercritical or uniform flow.SGS output is the product of the head drop over the dam and the flow ratethrough the dam. Therefore if the flow in Seton Canal is reduced due to hydroki-netic turbines, generation will be lost at SGS. If hydrokinetic turbines are to gen-erate electricity without diminishing SGS output, there must be uniform flow up-stream of the installation (condition b). The existence of uniform flow only guar-antees that an infinitely small amount of hydrokinetic generation could be installedwithout effect. Large amounts could overwhelm the control of both uniform andsupercritical flow.4.2 Model DevelopmentA one-dimensional steady state, finite element model was written in Matlab inorder to determine the effects of hydrokinetic turbines, The code may be found inAppendix D.The length of Seton Canal is discretized into 1m cells and is solved moving103Water Level:Before turbine installationAfter turbine installationHydrokineticTurbinesUniform, Subcritical FlowHydraulic JumpWater Level:Before turbine installationAfter turbine installationHydrokineticTurbinesSupercritical FlowSubcritical FlowWater Level:Before turbine installationAfter turbine installationAfter reduced flowHydrokineticTurbinesSubcritical FlowDiscontinuityReservoirControlledWater Level(c) Backwatering reaches controlled water level.(b) Backwatering into uniform, subcritical flow.(a) Backwatering into supercritical flow.Figure 4.2: Effect on surface profiles after installing hydrokinetic turbines104upstream, applying conservation of energy. See Figure 4.3.h1mFigure 4.3: Trapezoidal cells showing flow cross section (used to calculatevelocity) and potential head differenceThe developed model simulates nominal flow, taking into account canal geom-etry and changes in elevation, calibrated to field data onsite. This nominal output isvalidated against the available data. Both model and field data are used to confirmsubcriticality of the flow throughout the canal.Once validated, the effects of hydrokinetic turbines are added. The backwater-ing of the turbines propagates upstream; this distance and magnitude of backwater-ing is a function of the capacity of the hydrokinetic devices. If this backwateringreaches the model boundary, the corresponding effect on SGS generation is deter-mined.Many parameters were taken from engineering and survey drawings of Setoncanal, however these drawings cannot be released publicly.1054.2.1 Selection of Model BoundariesIt is important that the selection of model boundaries accurately represents theoperational constraints of the Seton System:• The downstream boundary is selected as the head-pond at station (the lastpoint of elevation survey data). Given that changes in Seton Canal are un-likely to affect the water level in the Fraser River at the SGS outlet, the SGSproduction is a function simply of the water surface elevation at the head-pond and the flow rate conveyed by the canal.• The upstream boundary is selected as the entrance to the canal at SetonReservoir. This boundary is selected for modeling efficiency as the sys-tem becomes significantly more complex if the reservoir level is allowed tochange. Such variation would have impacts on the outflow into the SetonRiver and land use on the shores of the reservoir. Furthermore, SGS genera-tion would be maximized if the reservoir is kept at its highest feasible levels;therefore that this is assumed to be the nominal state, thus hydrokinetic gen-eration should not result in an increase in reservoir level.This selection of boundaries provide a relatively simple model domain withrelatively few physical relationships needed to describe the behaviour of the flow.4.2.2 Theory and AssumptionsThe developed model is a numerical simulation that steps through the length of thecanal, solving for conservation of energy between each cell.Seton Canal, at 3.5 km and only 30 m wide, 4-6 m deep is modelled as a one-dimensional system as vertical and transverse motions are not considered to be106significant and are not of interest in this context (only longitudinal velocity is im-portant). This does not imply that the flow is one-dimensional, simply that effectssuch as velocity profiles are accounted for in the formulation of various parame-ters such as frictional resistance and are not simulated directly. This dramaticallyincreases the speed of the model, allowing it to be run over a large range of param-eters.The convention of mechanical energy in a flow being referred to as variousforms of head will be used in the following discussion.Energy BalanceThe energy in a fluid has a free exchange between elevation, pressure and velocityhead. As an open channel, the pressure at the free surface is always atmospheric(zero-gauge) and may be neglected. As a result, energy is converted between ele-vation and velocity head as (4.1)E =v22g+ y+ z (4.1)where v is the fluid velocity, y is the fluid depth, z is the elevation of the channelbottom and g is gravitational acceleration. Complicating this balance for flow inopen channels is that velocity and depth are interdependent. For a fixed volumetricflow rate, Q, velocity is given as (4.2)v =QA=Qy∫0b(y)dy(4.2)107where A is the cross sectional area of the flow and b(y) is a function that describesthe breadth of the channel. In calculation of a particular cell, the model iterates toensure the relationships between depth, velocity, flow and channel geometry are allsimultaneously satisfied.Energy in a model cell is calculated as the total-head in the adjacent cell down-stream, accounting for frictional losses and hydrokinetic generation as (4.3)Ei+1−h f −ht = Ei (4.3)where i is the cell index, h f is head lost due to friction and ht is head extractedby hydrokinetic turbines. The formulation of h f and ht will be discussed in thesections below.Frictional LoosesThe parameterisation of friction is quite important and a number of relationshipsexist. Friction is a function of the shear stress developed by the flow along a solidboundary; this is an interplay between the area of the wall and the velocity gradientat it.The formulations for frictional head loss are chosen based on the general char-acteristics of the flow. This model uses the Darcy - Weisbach relation:h f = f · LRh ·v22g(4.4)where f is the friction coefficient, a property of the bed roughness, L is the longitu-dinal length of the cell and Rh is the hydraulic radius of the flow at that cell. Theseparameters are calculated at each step before the change in energy is applied.108Hydrokinetic GenerationHydrokinetic turbines are simulated in the model as a two-step process where thepower generated is calculated using a turbine coefficient but the effect on the flowis determined using the head loss over the turbine.Power extracted by the hydrokinetic turbines, Pe, is calculated as (4.5)Pe = kt · 12 ·ρ · v3 ·At (4.5)where At is the swept area of the turbine and kt is a turbine resistance coefficient,assumed to be 0.25, as tested using TD turbines. (Bowman, 2010) This representsthe best-case generation as it assumes that the turbines always operate at their mostefficient point. Using Pt , ht is calculated using the following relationship (4.6):Pe = 0.5 ·g ·ρ ·ht ·Qt (4.6)The assumed factor of 0.5 is applied to account for the efficiency of the turbine(if the data is available, it may be calculated as the quotient of the coefficient of per-formance on coefficient of thrust). There is a question, however, what proportionof the flow passes through the turbines and what fraction is diverted around them.An analysis by Garrett and Cummins (2004) (in relation to free-stream turbines intidal flows) has demonstrated that the greatest potential extraction exists when theturbines occupy the entire cross-sectional area of the flow, thusQtQ= 1, a conditionthat will be used in this modelling.109Changes in SGS GenerationAs was discussed in Section 2.3.4, the propagation of backwatering to Seton Reser-voir will have negative impacts on SGS generation. In this event, SGS could adjustoperation in one of two ways:• Continue to operate with existing head but a reduced flow rate (loweredflow), or• Reduce generator load such that the head-pond water level drops; this reestab-lishes the pressure gradient along the canal and maintains the nominal flowrate. (lowered head)For the sake of computational efficiency, the lowered head approach will beused unless otherwise specified.4.2.3 Integration of Seton Canal ParametersThe specific values for channel geometry, flow and power were found from a vari-ety of sources, with varying degrees of precision, certainty and universality.Channel Cross Sectional AreaThe function for cross sectional area of the flow originates from values taken offthe construction drawings of the canal (BC Electric, 1953, 1954a,b). From thesedrawings, a number of canal shapes, widths and maximum-depths were found. Forthe purpose of the model, widths and shapes were linearly interpolated betweenknown points. The maximum allowable depth of flow (available freeboard) alongthe length of the channel is relatively uncertain and was not available from theprovided drawings. As such, this model does not perform a check to determine if110the flow overtops the bank.Flow DepthBottom elevation was taken from a 2010 survey of the canal (BC Hydro Instru-mentation Services, 2010). Depths of flow along the channel were determined by asurvey of the site using an Acoustic Doppler Current Profiler (ADCP). In additionto providing velocity data, the ADCP returns the depth measured by each of its fourbeams. Due to the narrow width and sloping sides of the channel, the ADCP wasoriented such that the horizontal projections of two beams were always parallel tothe flow; only these beams were used in depth determination.The ADCP was towed up the canal during the course of normal operation toproduce a series of depths along the canal past stations 1+145 through 4+150. Inorder to produce a profile along the length of the channel, a handheld GPS unitwas fitted to the ADCP to produce location mapped depth points. These werethen filtered and interpolated to produce a contiguous profile. See Figure 4.4. Theflow through the canal is thought to have changed slightly during the survey. Theeffects of this are not observable in the data, but this variability is considered in thesensitivity analysis. See Section 4.3.2Flow RateIn addition to the longitudinal survey for depth, transverse surveys were made todetermine the flow rate in the canal, simply as the integration of velocity over theflow area. This gave a range of flow rates between 109 m3/s and 112 m3/s. Op-erationally, it is thought that the maximum flow rate is approximately 100 m3/s,(Percel, 2010) as such, this uncertainty will be investigated in the sensitivity anal-111ysis of this parameter.Figure 4.4: Locations of ADCP Depth Survey Points. Image adapted fromGoogle Earth.4.2.4 Calibration and ValidationThere is a vast array of possible errors when developing and using numerical mod-els, including engine (mathematical) errors, incorrect input parameters, misinter-pretation of results, miscalibration, and numerical artefacts. It is immensely im-portant to ensure that model results match with reality.As was discussed at the beginning of this section, a simple model allows rel-atively simple calibration. Beyond the geometric properties of the channel, twovalues are required to predict the behaviour of the flow: the downstream depth andthe friction coefficient, f. These two parameters were optimized to reduce the de-viation between model and observed output with the results shown in Figure 4.5.The root-mean-squared error is less than 1% of depth, which indicates a good rep-resentation of the flow. The majority of this small error may be attributable to noisein the observation data and three-dimensional effects caused by sudden changes inbottom elevation. The sharp spikes in the depth data are the result of bottom eleva-tion change relative to datum, not changes in the free-surface elevation. Supporting112this conclusion, a correlation of model-output to observed-output (Figure 4.6) doesnot show any trends in error.0 500 1000 1500 2000 2500 3000 3500 400044.555.566.5Distance Downstream from Dam (m)Depth (m)  Model OutputObserved DepthFigure 4.5: Model and Observed Depth for 100 m3/s flow4.3 Results and DiscussionThe model is run with the parameters discussed above, a flow of 100 m3/s anda depth profile as shown in Figure 4.7. The depth of the flow increases nearlymonotonically moving downstream, indicating non-uniform flow. Froude numbersare between 0.12 and 0.2 over the length of the cannal. As such, it is not possible toinstall any amount of hydrokinetic generation that will have zero negative impacton SGS production, as backwatering will always propagate to some extent to Seton1135.2 5.3 5.4 5.5 5.6 5.7 5.8 5.9 6 6.15.25.45.65.866.2Model Depth (m)Observed Depth (m)  Correlation Point1:1Figure 4.6: Correlation of depth from model and observed depthReservoir.The area swept by turbines is incrementally increased, resulting in hydroki-netic generation ranging from zero to 64 kW. Turbines are initially placed at 3300m from the dam crest. This placement ensures the maximum available distance toabsorb backwatering from the turbine. Once the swept area of the turbine becameequal to the cross-sectional area of the flow, multiple turbines are added with thebalance of swept area placed at 3100, 2900 and 2700 from the dam crest, in se-quence. As illustration, the water surface profile of 64 kW of generation is given1140 500 1000 1500 2000 2500 3000 3500 4000228230232234236238Distance Downstream from Dam (m)Elevation Above GSC Datum (m)  Water SurfaceCanal BottomFigure 4.7: Water surface elevation at 100 m3/s flow without hydrokineticgenerationin Figure 4.8, showing a moderate increase in water level at Seton Reservoir.The increase in water level is a function of the magnitude of hydrokinetic gen-eration, as shown in Figure 4.9 a and b. As can be seen, the change in water levelis reduced at distances away from the turbines but does not recover completely be-fore Seton Reservoir. The relatively level slope of this curve is illustrative of SetonCanal’s efficiency at conveying its flow.4.3.1 Effects on Seton GenerationThe magnitude of this loss of production is now examined and compared againstthe gain in production from hydrokinetic turbines. For every simulated condition,there is a loss of generation at SGS, as shown in Figure 4.10. The reduction ingeneration from SGS exceeds hydrokinetic generation in all cases except wherehydrokinetic generation is less than 2 watts. As such, it is not possible to aug-1150 500 1000 1500 2000 2500 3000 3500 4000228230232234236238Distance Downstream from Dam (m)Elevation Above GSC Datum (m)  Water Surface with 358m2Swept Turbine Area (64kW)Nominal Water SurfaceCanal BottomFigure 4.8: Water surface elevations for nominal and 64 kW of hydrokineticgeneration.ment net production at Seton Canal with a free-stream hydrokinetic turbine at anyreasonable level of electricity production.This calculation of SGS generation reduction is using the Lowered Head ap-proach as discussed in Section 4.2.2. For comparison, the Lowered Flow approachis modelled, considering a 200 m2 turbine (35.4 kW). To eliminate the discontinu-ity at the entrance to the canal, the flow rate must be reduced from 100 m3/s to 90.9m3/s resulting in a 9.1% drop in SGS output, which equates to roughly 4.2 MW. SeeFigure 4.11. The reduction in flow also decreases the hydrokinetic output from 35kW to 27 kW. As such, if turbines are installed, it will likely be beneficial to adjustSGS operation such that the head pond level is reduced and flow is maintained.11610210310410510−410−310−210−1100Increase in Water Levelat Seton Canal Inlet (m)(b) Power Generated from Hydrokinetic Turbines (W)0 500 1000 1500 2000 2500 3000 3500 400000.020.040.060.080.1(a) Distance Downstream from Dam (m)Change in Depthfrom Nominal (m)50m2 (8.6kW)100m2 (17.4kW)150m2 (26.4kW)200m2 (35.4kW)250m2 (44.3kW)Figure 4.9: (a) Variation from natural depth as a function of distance for vary-ing levels of hydrokinetic generation. (b) Change in water level fromnatural as a function of hydrokinetic generation4.3.2 Sensitivity AnalysisA sensitivity analysis was performed on the two measured variables (flow rate anddepth) and the assumed parameter from Section 4.2.2, henceforth referred to as theturbine drag ratio. See Figure 4.12. Each scenario is compared to the nominalcondition of 100 m3/s, observed depths from the 2011 depth survey and a turbinedrag ratio of 0.5. The independent parameter under consideration is the loss of11710−310−210−110010110210310−1100101102103104105106Turbine Swept Area (m2)Power (W)  Loss of Generation at SGSProduction by Hydokinetic TurbinesFigure 4.10: Loss of SGS generation and gain in hydrokinetic generation.generation at SGS for 25 kW and 50 kW hydrokinetic production.Flow rates of 90, 110 and 120 m3/s were considered and the magnitude of SGSloss is shown to be relatively insensitive to flow. It is important to note that the lossis relative to zero-hydrokinetic baseline in that scenario as SGS generation doesvary with flow.Accounting for measurement errors, two scenarios considered±10% deviationin flow depth. While the changes in output were quite sensitive to depth, thisparameter is least subject to this magnitude of error as the ADCP is unlikely toproduce errors of the corresponding magnitude.1180 500 1000 1500 2000 2500 3000 3500 4000−0.0200.020.040.060.08Distance Downstream from Dam (m)Deviation from nominal depth (m)  100m3/s90.9m3/sFigure 4.11: Deviation from nominal depth with unadjusted and adjustedflow rate for a turbine of 200 m2 swept areaTurbine drag ratio is an assumed parameter (0.5 nominal) and relates the pro-duction to the drag on the turbine and would be equal to 1 for a perfectly efficient(unrealistic) turbine. Unsurprisingly, the losses at SGS are highly sensitive to thisparameter and it should be minimized in the design process for a particular turbineand associated structure.From this analysis, it is unlikely that variation of input parameters would pro-duce a scenario where the installation of hydrokinetic turbines would result in a netpositive change in power production.4.4 Conclusions and RecommendationsHydrokinetic devices offer opportunity to produce electricity by extracting the ki-netic energy from moving water, however the installation of these turbines willhave an impact on the flow in which they are installed. If turbines are installed1190.8 0.85 0.9 0.95 1 1.05 1.1 1.15 1.2 1.2530405060708090100110120Variation from NominalLoss of SGS Generation (kW)FlowTurbine Drag RatioMeasured Depth50kW Hydrokinetic25kW HydrokineticFigure 4.12: Sensitivity of flow, depth and turbine drag ratio to loss of gener-ation at Seton.downstream of hydraulic controls, or downstream of uniform flow where a wa-ter level increase is not problematic, the installation of turbines would not causechanges to the bulk flow. These conditions typically occur in high-slope or highbottom friction environments.In the case of Seton Canal, the flow is subcritical and the installation of a tur-bine will cause backwatering, the impact of which will propagate upstream to SetonReservoir. Certain flow features such as uniform flow, supercritical flow or weirsmay serve to halt the propagation, however none of these features exist at Seton.120The results of simulation show that an installation of hydrokinetic turbines willhave a detrimental effect on generation at the SGS. Furthermore, for all reasonablelevels of generation, the loss of SGS production are larger than the gains from thehydrokinetic generation.A sensitivity analysis performed on measured and assumed parameters demon-strates that reasonable variation from the nominal scenario does not significantlychange the negative impact that hydrokinetic turbines would have on SGS produc-tion.For the purposes of maximizing power output, hydrokinetic turbines are not afeasible proposition at this location.121Chapter 5Conclusions andRecommendationsSociety needs to produce more clean, renewable electricity, but there is debateabout where this energy will come from. While large hydro, solar and wind powerhave lead the development of renewable energy, significant political, academic,and business effort has been invested in developing tidal and hydrokinetic energydevices. These are most often free-stream turbines.Free-stream turbines by definition do not have a head difference imposed onthem by a penstock, they are only exposed to a moving fluid. There is a commonmisunderstanding of the operation of free-stream turbines however, that despite nothaving a head difference applied, a head difference across the turbine occurs andis developed by the turbine itself. This is a necessary condition for the turbine toextract any amount of power.The effects of this head difference can propagate away from the turbine. In tidal122flows, this acts to reduce the bulk flow rate. In hydrokinetic flows, the head changecan impact other users of the river. Unlike wind and solar (which presently convertonly small fractions of the total solar or wind resource) free-stream turbines arecapable of extracting a significant proportion of total energy resource under regularconditions. As a result, installation of free-stream turbines can dramatically alterthe flow, with repercussions on power generation as discussed below.5.1 Tidal PowerMany resource inventories for tidal power have simply summed the kinetic energyflux across a number of cross sections of flow in channels. This methodology isindefensible as it neglects to include the effect of flow retardation that the turbineswill impose. An analytical model was developed that shows the maximum amountof power that could be extracted for a tidal channel is 38% of the natural fluidpower and computed as:ηmax =PmaxρgQ0∆Htide= 0.38 (5.1)This 38% limit is the extractible energy from the flow. The generation of usableelectricity has some amount of inefficiency so the maximum amount of electricityproduction is somewhat less, dependent on the turbine and generation equipmentbeing used.It is unlikely however that the maximum power will ever be extracted from aparticular tidal channel. To reach this maximum extraction, the tidal flow must bereduced to 58% of natural. It is unlikely that such an impediment of tidal flowwould be permitted by government.123The production of electricity from a given turbine is a function of its efficiency,its swept area and the velocity of water it intercepts. Flow reduction due to turbineresistance reduces the overall flow in the channel, not simply in the wake of thatturbine. This flow reduction will reduce the amount of production by any turbinealready installed in that channel and can increase the cost of energy across an arrayof turbines.The development of tidal energy technology is driven primarily by the marketand companies developing turbines ultimately hope to make profit selling and in-stalling these devices. At present, the initial installations of tidal turbines do notproduce electricity that is competitive with the market rates. There is the expec-tation from government and industry however, that the levelized cost of electricitywill drop as the technology becomes more mature.This assumption of LCOE reductions has undergone some preliminary testsand found that there is uncertainty of its long-term validity. The reduction in tech-nology costs as a function of installed capacity is a known model using the his-tory of other industries and this reduction in costs acts to reduce the LCOE of tidalpower at small levels of deployment. Flow reduction that occurs through the instal-lation of turbines sharply decreases the amount of energy available to be generatedand will generally increase the cost of the LCOE. The interplay between these twofactors will likely cause the LCOE of tidal energy to initially drop (due to techno-logical learning) and then increase ever after (due to flow reduction). Dependingon the size of the total tidal energy resource and the future market price for electric-ity, it is not clear if tidal power will ever become competitive on the open energymarket. There may be opportunity in higher-priced niche markets, however.The basis of the above conclusions is an analytical model of a tidal channel.124This analytical model was subjected to validation with a physical model. The ex-perimental apparatus was constructed using a number of motion control elementsand sensors, retrofitted into an existing water flume. This type of physical modeldoes not appear to have been described previously in the literature. The flume wasinstrumented and a number of experimental runs were conducted. The results ofthese experimental runs indicate good agreement between the analytical model andthe physical system.5.2 Hydrokinetic PowerHydrokinetic power is associated with tidal power due to the similarities of the tur-bine technology as both utilize underwater free-stream turbines. The environmentsand the physical characteristics that determine the respective resource are, however,dramatically different. Unlike tidal flows, the flow rate in rivers is essentially fixed.The installation of a resistive turbine element will not reduce the steady-state flowrate, but will change the water levels in the river, either upstream or downstreamdepending on the flow criticality.The effect of the turbine diminishes as one moves away from the device due tofrictional dissipation; hydraulic controls can also minimize the far-field impacts ofhydrokinetic generation. As a result, if hydrokinetic generation is installed close tolocations that are sensitive to changes in water level such as canals or hydroelectricdams, the detrimental effect of the hydrokinetic turbines should be assessed.The installation of hydrokinetic turbines in Seton Canal was analysed using apurpose-written 1D numerical simulation. The proposed hydrokinetic installationwould be in close proximity to an existing hydroelectric generation facility. Theresults clearly demonstrated that for this location the installation of turbines would125produce a total net reduction in energy output when the impact on the hydroelectricfacility was considered.5.3 Future WorkThere is an opportunity for future research on the topic of resource assessmentof tidal and hydrokinetic power. Of particular interest is refinement of the LCOEprojections. This work necessitates a more complete global model of tidal energyconsidering flow reduction. The primary challenge to this effort is the availabil-ity of suitable environmental data, particularly the ∆Htide across channels and itsrelationship with flow.The physical modelling of tidal flows can be extended through additional ex-perimental runs and expanding the range of water levels tested. This requires mod-ification of some of the physical components such as the bottom friction element.The experimental runs tested static differences in head across the channel. Furtherexperimentation investigating flows with changing head differences. Such workrequires upgrades to the apparatus, but would likely not be possible without signif-icant increases in laboratory pumping capacity.126BibliographyAbu-Khader, M. M. (2009). Recent advances in nuclear power: A review. Progressin Nuclear Energy. 10Amaral, S., Perkins, N., Giza, D., and McMahon, B. (2011). Evaluation of fishinjury and mortality associated with hydrokinetic turbines. Alden Research Lab-oratory Inc., Holden, Massachusetts, Electric Power Research Institute (EPRI).18Arrhenius, S. (1896). On the influence of carbonic acid in the air upon the temper-ature of the ground. Philosophical Magazine Series 5, 41(251):237–276. 4Atwater, J. (2008). Limitations on tidal-in-stream power generation in a strait.Master’s thesis, University of British Columbia. 25, 136Atwater, J. F. and Lawrence, G. A. (2011). Regulatory, design and methodologicalimpacts in determining tidal-in-stream power resource potential. Energy Policy,39(3):1694 – 1698. iv, 1, 24, 27Bae, Y. H., Kim, K. O., and Choi, B. H. (2010). Lake Sihwa Tidal Power PlantProject. Ocean Engineering, 37(5-6):454 – 463. 29127BBC (2010). Severn barrage tidal energy scheme scrapped by huhne. 15BC Electric (1953). Seton creek development cayoosh creek crossing: Plan andsections (u-36260). 110BC Electric (1954a). Cayoosh creek crossing flume concrete outline (u-36263).110BC Electric (1954b). Seton creek development canal typical lining details (u-36226). 110BC Hydro Instrumentation Services (2010). Seton canal level surveys 3+350 to4+437: 2004 and 2010 vertical displacement profile (13953-brg-son). 111Berkowitz, R. (2010). Britain Tidal Power Plan Put on Hold in Statement on En-ergy Policy. Science Magazine, ScienceInsider. 29Bowman, C. (2010). Powetech labs test results. Personal Communication. 109Boyce, J. R. (2013). Prediction and inference in the hubbert-deffeyes peak oilmodel. The Energy Journal. 3Brooks, A., Lu, E., Reicher, D., Spirakis, C., and Weihl, B. (2010). Demanddispatch. Power and Energy Magazine, IEEE, 8(3):20–29. 11Bryden, I., Grinsted, T., and Melville, G. (2004). Assessing the potential of a sim-ple tidal channel to deliver useful energy. Applied Ocean Research, 26(5):198 –204. xiii, 25, 136, 137, 138Clark, N. A. (2006). Tidal Barrages and Birds. Ibis, 148:152–157. 24128Draper, S., Adcock, T. A. A., Borthwick, A. G. L., and Houlsby, G. T. (2014).Estimate of the tidal stream power resource of the Pentland Firth. RenewableEnergy, 63(c):650–657. 25Ernst & Young LLP (2010). Cost of financial support for wave, tidal range andtidal stream operations. 46Federal Energy Regulatory Comission (2015). Hydrokinetic projects. 99Fountain, H. (2013). How to charge millions of electric cars? not all at once. 8Garrett, C. and Cummins, P. (2004). Generating Power from Tidal Currents. Jour-nal of Waterway, Port, Coastal and Ocean Engineering, 130(3):114 – 118. 25,32, 109, 139, 140, 141, 142, 144Garrett, C. and Cummins, P. (2005). The power potential of tidal currents in chan-nels. In Proceedings of the Royal Society A: Mathematical, Physical and Engi-neering Sciences, volume 461, pages 2563–2572, London. The Royal Society.xiii, 24, 25, 30, 35, 139, 145, 146, 147Gill, A. B. (2005). Offshore renewable energy: Ecological implications of generat-ing electricity in the coastal zone. Journal of Applied Ecology, 42(4):605 – 615.15, 16, 24, 29Goldemberg, J., Baker, J. W., Ba-N’Daw, S., Khatib, H., Popescu, A., Viray, F. L.,Anderson, D., Holdren, J. P., Jefferson, M., Jochem, E., Nakicenovic´, Reddy,A. K., Rogner, H.-H., Smith, K. R., Turkenburg, W. C., and Williams, R. H.(2000). World Energy Assessment: Energy and the challenge of sustainability.United Nations Development Programme. 28129Hagerman, G. and Bedard, R. (2006). Massachusetts Tidal In-Stream Energy Con-version (TISEC): Survey and Characterization of Potential Project Sites. ElectricPower Research Institute, TP-003-MA Rev 1 Massachusetts Site Survey Report.29Hagerman, G., Bedard, R., and Polagye, B. (2005). Guidelines for PreliminaryEstimation of Power Production by Tidal In Stream (Current) Energy ConversionDevices. Electric Power Research Institute, TP-001-NA Rev 3 Guidelines forPreliminary Estimation of Power Production. 1, 24, 38, 91Hagerman, G., Fader, G., and Bedard, R. (2006). New Brunswick Tidal In-StreamEnergy Conversion (TISEC): Survey and Characterization of Potential ProjectSites. Electric Power Research Institute, TP-003-NB Rev 1 New BrunswickSite Survey Report. 29, 42Hammons, T. J. (1993). Tidal power. Proceedings of the IEEE, 81(3):419 – 433.29Hardisty, J. (2009). The Analysis of Tidal Power. John Wiley and Sons. 12Henderson, B. D. (1984). The application and misapplication of the experiencecurve. Journal of Business Strategy. 46Holt, J. A. and Eaton, D. J. (2008). Assessment of the energy potential and esti-mated costs of wind energy in british columbia. Technical report, Garrad HassanCanada for BC Hydro. 12International Energy Agency (2000). Experience curves for energy technologypolicy. 46130International Energy Agency (2012). Technology roadmap: High-efficiency, low-emissions coal-fired power generation. 9International Energy Agency (2014). World energy investment outlook. 3Khan, M., Bhuyan, G., Iqbal, M., and Quaicoe, J. (2009). Hydrokinetic en-ergy conversion systems and assessment of horizontal and vertical axis turbinesfor river and tidal applications: A technology status review. Applied Energy,86(10):1823–1835. 16, 23Kim, K.-P., Ahmed, M. R., and Lee, Y.-H. (2012). Efficiency improvement ofa tidal current turbine utilizing a larger area of channel. Renewable Energy,48(0):557–564. 23Kummu, M., de Moel, H., Ward, P. J., and Varis, O. (2011). How close do we liveto water? a global analysis of population. PLoS ONE. 18Lutz, B. D., Bernhardt, E. S., and Schlesinger, W. H. (2013). The environmentalprice tag on a ton of mountaintop removal coal. PloS one, 8(9):e73203. 4Marine Current Turbines (2013). Seagen Technology; accessed August 7, 2013.http://www.marineturbines.com/SeaGen-Technology. 45Michal, R. (2001). Fifty years ago in december: Atomic reactor ebr-i producedfirst electricity. Nuclear News, 44(12):28–29. 9Myers, L. E. and Bahaj, A. S. (2010). Experimental analysis of the flow fieldaround horizontal axis tidal turbines by use of scale mesh disk rotor simulators.Ocean Engineering, 37(2–3):218–227. 23131Nuclear Energy Agency (2012). Country profile: France - summary figures for2012. 9Nyamdash, B., Denny, E., and O’Malley, M. (2010). The viability of balancingwind generation with large scale energy storage. Energy Policy, 38(11):7200 –7208. 11Percel, D. (2010). BC Hydro SGS Operation Parameters. Personal Communica-tion. 101, 111Polagye, B. and Bedard, R. (2006). Tidal In-Stream Energy Resource AssessmentFor Southeast Alaska. Electric Power Research Institute, TP-006-AK AlaskaTidal Power System Level Design. 1, 24, 29Previsic, M., Polagye, B., and Bedard, R. (2006). System Level Design, Perfor-mance, Cost and Economic Assessment – San Francisco Tidal In-Stream PowerPlant. Electric Power Research Institute, TP-006-CA California Tidal PowerSystem Level Design. 29Ramos, V., Carballo, R., A´lvarez, M., Sanchez, M., and Iglesias, G. (2013). As-sessment of the impacts of tidal stream energy through high-resolution numeri-cal modeling. Energy, 61(c):541–554. 25Rourke, F. O., Boyle, F., and Reynolds, A. (2010). Marine current energy devices:Current status and possible future applications in ireland. Renewable and Sus-tainable Energy Reviews. 23Royal Academy of Engineering (2004). The Cost of Generating Electricity. 41132Shafiee, S. and Topal, E. (2009). When will fossil fuel reserves be diminished?Energy Policy, 37(1):181 – 189. 3Shields, M. A., Dillon, L. J., Woolf, D. K., and Ford, A. T. (2009). Strategicpriorities for assessing ecological impacts of marine renewable energy devicesin the pentland firth (scotland, uk). Marine Policy, 33(4):635–642. 15Short, W., Packey, D. J., and Holt, T. (1995). A Manual for the Economic Eval-uation of Energy Efficiency and Renewable Energy Technologies. National Re-newable Energy Laboratory, U.S. Department of Energy. 45Sims, R., Schock, R., Adegbululgbe, A., Fenhann, J., Konstantinaviciute, I.,Moomaw, W., Nimir, H., Schlamadinger, B., Torres-Martı´nez, J., Turner, C.,Uchiyama, Y., Vuori, S., Wamukonya, N., and Zhang, X. (2007). Energy supply.In Climate Change 2007: Mitigation. Contribution of Working Group III to theFourth Assessment Report of the Intergovernmental Panel on Climate Change.Cambridge University Press, Cambridge, United Kingdom and New York, NY,USA. 5, 6Solomon, S., Qin, D., Manning, M., Chen, Z., Marquis, M., Averyt, K., Tignor,M., and Miller, H., editors (2007). Contribution of Working Group I to theFourth Assessment Report of the Intergovernmental Panel on Climate Change.Cambridge University Press. xi, 6Stern, D. I. and Cleveland, C. J. (2004). Energy and economic growth. Rensselaerworking papers in economics, Rensselaer Polytechnic Institute, Department ofEconomics. 2133Tarbotton, M. and Larson, M. (May 2006). Canada Ocean Energy Atlas (Phase 1)Potential Tidal Current Energy Resources. Report for the Canadian HydraulicsCentre, Ottawa. 1, 15, 24, 29, 30, 37, 48, 148Tidmarsh, W. G. (1983). Assessing the environmental impact of the annapolis tidalpower project. Water Science and Technology, 16(1-2):307 – 317. 29Triton Consultants Ltd. (2002). BC Green Energy Study Phase 2: Tidal CurrentEnergy. Report for British Columbia Hydro and Power Authority. 29Trussart, S., Messier, D., Roquet, V., and Aki, S. (2002). Hydropower projects:a review of most effective mitigation measures. Energy Policy, 30(14):1251 –1259. 10US Department of Commerce, NOAA, Earth System Research Laboratory, andESRL Web Team (2014). ESRL Global Monitoring Division - Global Green-house Gas Reference Network; accessed March 1, 2014. 5US Energy Information Administration. Renewable Energy Sources; accessedSeptember 27, 2013. 11US Energy Information Administration (2012). Regulators approve first commer-cial hydrokinetic projects in the United States; accessed September 20, 2014.http://www.eia.gov/todayinenergy/detail.cfm?id=8210. 100US Energy Information Administration (2013a). Assumptionsto the Annual Energy Outlook; accessed September 20, 2014.http://www.eia.gov/forecasts/archive/aeo13/. 46134US Energy Information Administration (2013b). International Energy Outlook2013 with Projections to 2040; accessed September 21, 2014. 8US Energy Information Administration (2014). International Energy Outlook 2014with Projections to 2040; accessed September 21, 2014. 7US Environmental Protection Agency. Radioactive Waste Disposal: An Environ-mental Perspective; accessed March 20, 2013. 10van den Wall Bake, J., Junginger, M., Faaij, A., Poot, T., and Walter, A. (2009).Explaining the experience curve: Cost reductions of brazilian ethanol from sug-arcane. Biomass and Bioenergy, 33(4):644 – 658. 46Vennell, R. (2012). The energetics of large tidal turbine arrays. Renewable Energy,48(0):210–219. 23Walters, R. A., Tarbotton, M. R., and Hiles, C. E. (2013). Estimation of tidal powerpotential. Renewable Energy, 51(0):255–262. 25World Bank (2014). World Development Indicators: Electricityproduction, sources, and access; accessed September 20, 2014.http://data.worldbank.org/indicator/EG.ELC.HYRO.ZS. 2, 9World Energy Council (2013). World Energy Resources: Marine Energy; accessedMarch 20, 2013. https://www.worldenergy.org/data/resources/resource/marine/.45135Appendix ALiterature Review of GeneralAssessment MethodologiesNote: This section is taken directly from Atwater (2008) to be used as reference.An analysis by Bryden et al. (2004) uses a momentum balance to predict thechanges in flow with small energy extractions. The authors consider a slug of fluidmoving through a channel. Assuming steady state, the pressure gradient producedby the tidal head is balanced by advection and a retarding force of natural frictionand the extraction devices such that:pressure gradient︷ ︸︸ ︷ρgUA∂h∂x=−advection︷ ︸︸ ︷U∂∂x(AU2ρ)−bottom friction︷ ︸︸ ︷UPerτ0 −power extraction︷︸︸︷PxA (A.1)where U is the fluid velocity, A is the cross sectional area, Per is the wettedperimeter, τ0 is the natural shear stress and Px is the power extraction per unitvolume (W/m3).136Applying conservation of mass for steady flow with no lateral inflow, such that∂Q∂x = 0 Equation (A.1) can be expressed as:(1− Q2h3b2g)∂h∂x=∂b∂xQ2gh2b3− 1ρgbhPerτ0− PxAρgQ (A.2)Bryden et al. (2004) then integrate Equation (A.2) with Px = 0 to determinethe surface profile for an arbitrary Q. Various solutions of (A.2) are iterated to de-termine the flow that corresponds to the tidal head under investigation. The powerextraction term is then increased and the relationship solved as before (Figure A.1).The author ends the analysis at this point; the primary goal of the paper was show-ing that the flow does indeed slow when additional retarding force is added (FigureA.2). The relevant figure in Bryden et al. (2004) does not show that the flow quicklydeviates from its initial linear trend at high levels of extraction. Using an arbitrarychannel geometry, the authors calculate that extracting 10% of the KE flux willreduce the flow by 3% and extracting 20% will reduce the flow by 6%. It is thensuggested that extraction of 10% of the KE flux could be considered a ‘rule ofthumb’ for future development of tidal power. Bryden et al. (2004) do not provideany rationale for the value of 10% extraction.The formulation of the momentum equation used in this paper does not lenditself to the estimation of total extractible power as Equation (A.2) relies on a Px torepresent the turbines.The magnitude of power extraction is a function of the flow and the turbinecharacteristics. In formulation used by Bryden et al. (2004) however, the powerextraction term is a control for the flow. In order to determine the maximum ex-traction, (a topic not addressed by the authors) one must slowly increment Px until13739.139.239.339.439.539.639.739.839.940.040.10 500 1000 1500 2000 2500 3000 3500 4000Distance (m)Depth (m)0.000.501.001.502.002.50Velocity (m/s)DepthVelocityLocation of ConversionDevicesFigure A.1: Channel flow properties with no artificial energy extraction.Recreated from Bryden et al. (2004)-50%-45%-40%-35%-30%-25%-20%-15%-10%-5%0%0% 5% 10% 15% 20% 25% 30% 35% 40%Influence of Energy Extraction on Current SpeedChange in flow speedShown in publicationNot Shown in publicationFigure A.2: Influence of proportional extraction on the mean flow speed.Recreated with modification from Bryden et al. (2004)the model fails. Furthermore, since the controlling parameter is dimensional, it isdifficult to apply the results to multiple scenarios. In essence, this parameter is notparticularly useful in the determination of maximum power extraction in the gen-eral case. However it can be used to determine the effects of small, fixed amountsof extraction.Despite the limitations of using the Bryden et al. (2004) approach for the es-138timation of total extractable power, a very important aspect is addressed in thatthe flow speed will be reduced with the installation of extraction devices. Further-more, the authors remind readers that the design of extraction devices is largelydependent on the flow velocities for the determination of optimum geometry. Tidalenergy developers must take into account the reduced velocities and not base theirdesigns on the natural conditions.More rigorous analysis that addresses total extractible power was conductedby Garrett and Cummins (2004) and Garrett and Cummins (2005). Garrett andCummins (2004) provides some key insights as to the behaviour of tidal flows intoand out of bays whileGarrett and Cummins (2005) addresses flows in straits.The first point presented by Garrett and Cummins (2004) is to remind the readerthat for a turbine to extract power it must experience a pressure difference and musthave flow passing through it; the absence of either of these conditions will result inzero power production. In order to perform work, a fluid requires both flow and adifferential pressure as shown by:W = F ·d→ dWdt=∫p~u · nˆdA (A.3)where W is work, F is force, d is distance, t is time, u is fluid velocity, p is pressureacting upon the surface under consideration and A and nˆ the area and normal vectorrespectively of that surface.Garrett and Cummins (2004) argue that maximum power is produced when theturbines span the entire cross section of the tidal channel. A general situation ispresented where a turbine is installed in a flow such that it covers some fraction oftotal cross sectional area, denoted by ε (Figure A.3).139u0 T p0 p2 p3p1u1 u2 u0Figure A.3: Schematic of flow through turbine or turbines, denoted by T, oc-cupying fraction ε of cross-sectional area of channel. Recreated fromGarrett and Cummins (2004)Assuming steady flow and a constant cross section, continuity then requiresεu1+(1− ε)u2 = u0 (A.4)and thusu2 =u0− εu11− ε (A.5)where the subscript 0 denotes the region upstream of the turbine, 1 denotes theregion through the turbine, 2 is the region to either side of the turbine and 3 is theregion downstream of the turbine.Garrett and Cummins (2004) assume that there is not a significant transversepressure gradient downstream of the turbine. This assumption is likely reasonablein the case of a turbine with a low blockage factor. Using Bernoulli’s theorem, thepressure drop across the turbines is then given by:p0− p2 = 12ρ(u22−u20) (A.6)140and thus power given byP =12ερu1(u22−u20) (A.7)Using Equation (A.7), power is compared to the situation where the turbinesoccupy the entire channel, such that:P0 = (p0− p3)u0 (A.8)This yields:PP0=x(1+ x−2εx)2x+ ε(1−3x) (A.9)where x = u1/u0. Given that if the turbine is to act as a turbine and not a pump,x must be be less than unity. Solving Equation (A.9) for 0 < x < 1 (Figure A.4)shows that P/P0 ≤ 1 for all values of ε . This implies that power from a set ofturbines occupying a fraction of the channel cross section will never be larger thanfrom a set of turbines occupying the entire channel cross section. As such, it canbe concluded that in determining the maximum power output, only the situation inwhich the entire channel area is impeded need be considered.Garrett and Cummins (2004) then focus on extractible power through a channelleading to a closed bay. In this scenario the difference in water level between theinside of the bay and the ocean drives a flow. This flow is resisted by a turbineforce per unit mass, F , creating a dynamical balance given by:g(ζ −acosωt)L=−F (A.10)where g is gravitational acceleration, a and ω are the amplitude and frequency oftidal forcing respectively, t is time, L is the length of the channel and ζ is the water14100.10.20.30.40.50.60.70.80.910 0.2 0.4 0.6 0.8 1xP/P0!=0.9!=0.5!=0.1Figure A.4: Solution of Equation (A.9) for ε = {0.1,0.5,0.9}. Recreated from Garrett and Cummins (2004).surface elevation inside the bay. Two models for friction are then applied, onelinear with velocity and one quadratic with velocity such that F = r1u or F = r2u|u|.Continuity is applied to the channel/bay as:Adζdt= Eu (A.11)where E is the cross sectional area of the channel and u is the velocity in thechannel. Equations (A.10) and (A.11) are non-dimensionalised and combined with142ζ = aζ ‘, t = t‘/ω and u = (Aωa/E)u‘ giving:dζ ‘dt‘=cos t‘−ζ ‘δ1(A.12)for the linear case where the non-dimensional resistance is:δ1 =r1LAωgE(A.13)and for the quadratic case:dζ ‘dt‘=|cos t‘−ζ ‘|1/2sgn(cos t‘−ζ ′)δ2(A.14)where the non-dimensional resistance is:δ2 =r2LA2ω2agE2(A.15)Solving the dynamical balance given in Equation (A.10) allows the calculationof power production using P = ρELFu givingP =(2δ11+δ 21)Pmax (A.16)for linear resistance andP = 4δ2∣∣∣∣(dζ ‘dt‘)∣∣∣∣Pmax (A.17)for quadratic.As can be seen in Figure A.5, power production rises to a maximum as turbineresistance is added. Past this peak, production falls with added turbine resistance14300.10.20.30.40.50.60.70.80.910 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5!1P/PmaxFigure A.5: Ratio P/Pmax as function of δ1 (linear resistance). Recreatedfrom Garrett and Cummins (2004)as the effect of the flow slowing becomes the dominant factor.An interesting additional conclusion is that the tidal range inside the bay is notnecessarily significantly affected even at the maximum extraction rate; the headacross the turbines is largely driven by the phase difference of the water levelsinside and outside the bay rather than solely by the differences in sinusoidal ampli-tudes. As such, flushing of the bay can be maintained.A specific note is made by Garrett and Cummins (2004) as to the use of ki-netic energy flux as a metric for available power, specifically stating: “Perhaps thesafest conclusion is just that 1/2ρu3 is unreliable as a power estimate in a situationinvolving flow into a closed basin”The conclusions of this article have some significant implications with regardsto resource assessment of tidal energy:144• The generation of power from a turbine requires both flow and differentialpressure.• Analysis may be simplified as maximum generation occurs when all of aflow passes through the turbines.• For a closed bay, generation does not necessarily significantly reduce tidalflushing due to a phase difference that develops between the inner and outerwater levels. It is important to note that this last conclusion can not be ap-plied to tidal straits (the focus of this thesis) as the head difference that drivesthe flow is dependent on ocean water levels and is not generally changed bythe flow through the channel.Garrett and Cummins (2005) address the issue of generation in tidal straits using asimilar approach. The momentum equation for this case is given by∂u∂ t+u∂u∂x+g∂ζ∂x=−F (A.18)where the pressure gradient along the channel is balanced by some element ofresistance, advection and acceleration. Garrett and Cummins (2005) then subjectsthis relationship to a variety of conditions.The initial test is to neglect natural friction and to balance the pressure gradientonly with acceleration and turbine resistance (the authors acknowledge that thissituation is unlikely to occur in reality.) As in Garrett and Cummins (2004) powergeneration is given by ρuAF , now integrated over the length of the channel and145assuming a linear turbine drag givingP =∫ L0ρFQdx = ρQ∫ L0Fdx (A.19)yieldingcdQdt−gacosωt =−λQ (A.20)where λ is a term associated with the number of the turbines in the flow.Assuming that the flow is uniform across each cross-section, Equation (A.20)can be solved, giving the maximum extractible power asP0 =14ρg2a2(cω)−1 =14ρcωQ20 =14ρgaQ0 (A.21)where a is the amplitude of the tidal head across the channel and Q0 being the nat-ural flow rate. While a unique solution may be reached, the fundamental assump-tion that bottom friction is negligible is unrealistic as natural friction is dominantin many situations.Garrett and Cummins (2005) continue with analysis neglecting bottom frictionbut investigating non-linear turbine resistance with flow. It was found throughnumerical solutions of the momentum equation that maximum extractible poweroccurs with linear friction. Various power extraction values are presented. SeeFigure A.6. It should be noted that resistance to flow is typically considered to bequadratic with velocity (n = 2)Assuming the drag is quadratic with flow, the authors add the effects of bottomfriction and exit losses from the channel. The friction term, λ , is separated intoa turbine resistance and a friction resistance; the ratio of power generated by the146000.20.40.60.81.01 2 3 4 5P/P0λ’n=0.25n=0.5n=1n=1.5n=2n=3Figure A.6: The scaled maximum power as a function of a parameter λ ′ ,representing the frictional drag associated with the turbines, for vari-ous values of n where the turbine drag is assumed proportional to thenth power of the current speed. Recreated from Garrett and Cummins(2005)turbines to that dissipated by friction is equal to the ratio between the separatedturbine friction parameter and the total friction parameter used in Equation (A.20).Considering situations where the natural friction coefficient is large, (thus fric-tion is dominant, not acceleration), the ratio of extractible power vs friction-freeextractible power approaches 0.86λ−0.5. Using this result, the extractible powerfrom a strait can be given by 0.38ρgQ1ζ . Where Q1 is the natural flow rate and ζis the head difference across the channel. This corresponds to a reduction in flowrate to 58% of natural.Garrett and Cummins (2005) make a thorough argument as to the maximum of147extractible power using a fundamental momentum balance. While the approach iscomplete, it has two weaknesses:• it is somewhat difficult to understand the procedure.• the model is not easily adaptable for variable turbine behaviour.Present assessments of tidal energy have been quite limited and have been in-tended to develop a rough estimate of the resource. As such, some studies havelacked the scope to investigate the complexities of tidal extraction and have electedto use the kinetic energy flux method due to its ease of application, a topic specif-ically addressed by Tarbotton and Larson (2006). The highly mathematical ap-proach used by Garrett and Cummins likely further inhibits adaptation of morecomplex assessment models in this context due to the time required to apply themethodology.In general, Garrett and Cummins’ arguments are generally complete for thecases they addressed, however they may be limited in application due to its highlymathematical nature and lack of easy adaptability.148Appendix BPython Code for Calculation ofTidal Power LCOEimpor t pandas as pdfrom numpy impor t *from f l u i d c o n s t s impor t *def OC wi th learn ing (C, L ) :C0 = L [ 'C0 ' ] # i n i t i a l cumulat ive i n s t a l l a t i o n sOC = 1.0#Per iod 1b = - log (1 -L [ ' l r 1 ' ] ) / log ( 2 . 0 )a = OC / L [ 'C0 ' ] * * ( - b )i f C <= 2**L [ ' p1 ' ] :OC = a * C * * ( - b )e lse :OC = a * ( ( 2 * * L [ ' p1 ' ] ) * * ( - b ) )#Per iod 2b = - log (1 -L [ ' l r 2 ' ] ) / log ( 2 . 0 )a = OC / ( ( 2 * * L [ ' p1 ' ] ) * * ( - b ) )i f C <= 2 * * ( L [ ' p2 ' ] ) :OC = a * C * * ( - b )e lse :149OC = a * ( 2 * * L [ ' p2 ' ] ) * * ( - b )#Per iod 3b = - log (1 -L [ ' l r 3 ' ] ) / log ( 2 . 0 )a = OC / ( ( 2 * * L [ ' p2 ' ] ) * * ( - b ) )OC = a * (C * * ( - b ) )r e t u r n OCdef LCOE(CapEx , OpEx ,AEP, L i fe , i ) :CRF = ( i * (1 + i ) * * L i f e ) / ( ( ( 1 + i ) * * L i f e ) -1)p r i n t CRFLCOE = (CapEx*CRF + OpEx) / AEPr e t u r n LCOELCOE(13750000*0.29 ,300000*.29 ,6000 ,10 , .15)OC = pd . DataFrame ( range (0 ,10000) )OC=OC. rename ( columns ={0: 'N ' })OC[ 'OC ' ] = NoneOC[ 'Cuml OC ' ] = NoneOCsum = 0f o r i i n OC. index :i f i > 0:i f i < 100E3 :OC. i x [ i , 'OC ' ] = OC wi th learn ing (OC. i x [ i , 'N ' ] , l ea rn i ng )OCsum = OCsum+OC. i x [ i , 'OC ' ]OC. i x [ i , 'Cuml OC ' ] = OCsumdef Power ( t , c , N, expanded data=False ) :# t : Turbine D ic t# c : Channel D i c t# N: Number o f t u rb i nes i n s t a l l e d# Returns the e l e c t r i c i t y ex t rac ted from the channel i n Wattsdebug = 0#Assign to more readable v a r i a b l e s#Channelf = c [ ' f ' ]u0 = c [ ' u0 ' ]A c = c [ ' Depth ' ] * c [ ' Breadth ' ]Rh c = 4* A c / (2 * c [ ' Depth ' ] + c [ ' Breadth ' ] )Q0 = u0* A cL c = c [ ' Length ' ]#TurbineD t = t [ ' Diameter ' ]A t = p i * ( D t * * 2 ) / 4cp = t [ ' cp ' ]c t = cp * t [ ' c t on cp ' ]150#Convert to k fk f = f * ( L c / Rh c ) / (2 * g * A c * * 2 )# Ca lcu la te the t o t a l head dropdH = k f * Q0**2i f debug >= 6 : p r i n t 'dH : ' ,dH, 'm '# Ca lcu la te k tk t = (N * c t * A t ) / (2 * g * A c * * 3 )i f debug >= 5 : p r i n t ' k f : ' , k f , ' k t : ' , k t ,# Ca lcu la te new f low ra teQ = (dH / ( k t + k f ) ) * *0 .5i f debug >= 5 : p r i n t 'Q: ' ,QV = Q/ A cPex t rac t = rho * g * k t * Q**3Pgen = Pex t rac t * ( cp / c t ) *0.556i f debug >= 2 : p r i n t 'V : ' ,V , ' Pgen : ' ,Pgenoutput = {' Pgen ' : Pgen ,' k t ' : k t ,'V ' :V ,'N ' :N,}i f expanded data :output [ 'Q ' ] = Qoutput [ 'dH ' ] = dHoutput [ ' k t o n k f ' ] = k t / k foutput [ ' eta ' ] =Pgen / ( rho * g * Q0 * dH)output [ 'Q ' ] =Qoutput [ ' Pex t rac t ' ] =Pex t rac toutput [ ' k f ' ] = k fr e t u r n output# <codecel l>f i n a n c i a l = {' i r a t e ' : 0 .15 ,' L i f e ' : 10 ,}t u r b i n e = { ' cp ' : 0 . 35 ,' c t on cp ' : 2 ,' Cost ' : 1E6 ,' OpExFraction ' : 0 .05 ,' Diameter ' : 20 ,' Breadth diameters ' : 2 ,' Stream diameters ' : 10151}channel = {' f ' : . 0 1 ,' u0 ' : 2 . 0 , #m/ s , t h i s w i l l be changed i n the loop' Depth ' : 50 .0 , #m' Breadth ' :400 , #m' Length ' :20000.0 , #m}l ea rn i ng = { ' p1 ' : 3 ,' p2 ' : 8 ,' l r 1 ' : 0 .2 ,' l r 2 ' : 0 .1 ,' l r 3 ' : 0 .01 ,'C0 ' : 1 ,'max ' : 1E6 ,}# <codecel l># loop through v e l o c i t i e sranges = {' l ow res ' :100E6 ,' h igh res ' :250E6 ,' res s tep ' : 100E6 ,' l ow ve l ' : 1 . 0 ,' h i g h v e l ' : 3 ,' v e l s t e p ' : 0 . 5 ,}# load the channel l i s t from Tarbot ton Datachanne l i n fo = pd . read exce l ( ' Channels . x l sx ' , ' l a r g e s t ' )c h a n n e l l i s t = [ ]f o r i i n range (0 , len ( channe l i n fo ) ) :c = {' f ' : . 002 ,' u0 ' : channe l i n fo . i x [ i , ' V e l o c i t y ' ] , #m/ s , t h i s w i l l be changed i n theloop' Depth ' : channe l i n fo . i x [ i , ' Depth (m) ' ] , #m' Breadth ' : channe l i n fo . i x [ i , ' Width (m) ' ] , #m' Length ' : channe l i n fo . i x [ i , ' Width (m) ' ] * 5 , #m' A t ' : 0 ,' u ' : channe l i n fo . i x [ i , ' V e l o c i t y ' ] ,'N ' : 0 ,' Power ' : 0}c h a n n e l l i s t . append ( c )#Change the channel leng th to vary the resource p o t e n t i a ldef c a l c p o t e n t i a l ( c h a n n e l l i s t ) :ke = [ ]152po ten t i a l power = [ ]f o r i , ch i n enumerate ( c h a n n e l l i s t ) :A c = ch [ ' Depth ' ] * ch [ ' Breadth ' ]Rh c = 4* A c / (2 * ch [ ' Depth ' ] + ch [ ' Breadth ' ] )k f = ch [ ' f ' ] * ( ch [ ' Length ' ] / Rh c ) / (2 * g * A c * * 2 )Q = A c * ch [ ' u0 ' ]dH = k f * Q**2po ten t i a l power . append (0.221 * rho * g * dH * Q)ke . append (0 .424*0 .5 * rho * A c * ch [ ' u0 ' ] * * 3 )r e t u r n po ten t ia l power , keo u t p u t l i s t = [ ] #This i s the l i s t where the output dataframes w i l l be s toredf o r desired AEP i n [1 ,10 ,100 ] : #Desired annual energy produc t ion i n TWh/ yr# p r i n t d i s i r e d p o t e n t i a l , 'MW, ( ' , desired AEP , ' TWh/ y r ) '#conver t to d i r e c t mean i n s t a l l e d power (W)d i s i r e d p o t e n t i a l = ( desired AEP *1E12/8760)# Ca lcu la te the maximum amount o f power a v a i l a b l e using 0.221* rho *g*dH*Q0# I t assumes the channel leng th i s = channel breadthpo ten t ia l power , ke = c a l c p o t e n t i a l ( c h a n n e l l i s t )#So we now change the channel leng th by a f a c t o r to get the r i g h t powerp o t e n t i a lm u l t i p l i e r = d i s i r e d p o t e n t i a l / np . ar ray ( po ten t i a l power ) . sum ( )f o r i , ch i n enumerate ( c h a n n e l l i s t ) :c h a n n e l l i s t [ i ] [ ' Length ' ] = c h a n n e l l i s t [ i ] [ ' Length ' ] * m u l t i p l i e r#Grab the new p o t e n t i a l to checkpo ten t ia l power , ke = c a l c p o t e n t i a l ( c h a n n e l l i s t )#Output the undis turbed channel data i n t o a dataframe to add tu rb i nes tod f = pd . DataFrame ( c h a n n e l l i s t )# In the loop , one t u r b i n e i s added at a t ime . This i s where the power outputgoes .t o ta l power = [ ]N=0 #Number o f i n s t a l l e d tu rb ined . ( Counter )# I t w i l l keep going u n t i l there are no channels w i th f low above 1.5m/ s . Thisi s a lower l i m i t on the p o t e n t i a l generat ion .p r i n t ' S t a r t i n g a d d i t i o n o f t u rb i nes a t ' , d i s i r e d p o t e n t i a l , 'MW ( ' , desired AEP, ' TWh/ y r ) 'wh i le d f [ ' u ' ] . max ( ) >1.5:# increment counterN=N+1# f i n d maximum v e l o c i t y o f the channels . This i s where the next t u r b i n ew i l l be i n s t a l l e d153i = d f . i x [ d f . u == df . u .max ( ) ] . index [ 0 ]#Put the channel w i th max v e l o c i t y i n t o a d i c t ob jec t#This i s here because of the prev ious vers ion o f code took a d i c t r a t h e rthan a dataframe as inpu t .#We conver t back and f o u r t h so I didn ' t have to r e w r i t e the Power ( ) code .th i s chan = d i c t ( d f . i x [ i ] )#Send the d i c t to the Power f u c t i o n . I t r e tu rns a d i c t w i th the newv e l o c i t y and generated power ( i n e l e c t r i c i t y )output = Power ( tu rb ine , th is chan , th i s chan [ 'N ' ]+1 )#Put the output data i n t o the d i c t i o n a r y t h a t can be shoved back i n t o thedataframe .th i s chan [ ' Power ' ] = output [ ' Pgen ' ]t h i s chan [ 'N ' ] = output [ 'N ' ]t h i s chan [ ' u ' ] = output [ 'V ' ]#Shove i t back i n t o the dataframe .d f . i x [ i ] = pd . Ser ies ( th i s chan )#Put the cumulat ive power i n t o a l i s t o f d i c t s so i t can e a s i l y beconverted i n t o a dataframe .tp = { ' Power ' : d f [ ' Power ' ] . sum ( ) ,'N ' :N,}t o t a l power . append ( tp )i f N%1000 == 0:p r i n t Np r i n t ' [DONE] ' , len ( t o ta l power ) , ' t u rb i nes i n s t a l l e d . '#Convert the l i s t o f d i c t s o f power output i n t o a dataframe f o r f u r t h e r workpdf=pd . DataFrame ( to ta l power )#Set the nominal nameplate i n MW. This i s taken from MCTturb ine namepla te = 1.2#Add columns to the dataframepdf [ ' t o t a l c a p e x ' ] = Nonepdf [ ' t o t a l o p e x ' ] = Nonepdf [ 'AEP ' ] = None #MWhpdf [ 'LCOE ' ] = None #$ /MWh#Opex i s given per MW from the Ernst and Young r e p o r topex per tu rb ine = 300E3* turb ine namepla te#Capex i s taken from MCT ' s number f o r Kyle Rheacapex per tu rb ine = 13.75E6#Loop through each increamenta l i n s t a l l a t i o n o f t u rb i nes . Ca lcu la te the LCOEat t h i s number o f t u rb i nesp r i n t ' S t a r t i n g c a l c u l a t i o n o f LCOE 'f o r i i n pdf . index :# I t s t ime consuming to c a l c u l a t e the f u l l amount o f cumulat ive overn igh tcosts and i t doesn ' t change much f o r la rge numbers .154#The OC i s ca l cu la ted f o r the f i r s t 10 ,000 tu rb ines , then l e v e l s out .i f pdf . i x [ i , 'N ' ] > OC. index [ - 1 ] :N = pdf . i x [ i , 'N ' ] #Number o f t u rb i nes i n t h i s loopabove max turb = N - OC. index [ - 1 ] #Number o f t u rb i nes above the maxd i r e c t l y computed# Ca lcu l t e the cucuml OC = OC. i x [OC. index [ - 1 ] , 'Cuml OC ' ] + (OC. i x [OC. index [ - 1 ] , 'OC ' ] *above max turb )pdf . i x [ i , ' t o t a l c a p e x ' ] = cuml OC* capex per tu rb inepdf . i x [ i , ' t o t a l o p e x ' ] = cuml OC* opex per tu rb inee lse :pdf . i x [ i , ' t o t a l c a p e x ' ] = OC. i x [ pdf . i x [ i , 'N ' ] , 'Cuml OC ' ] *capex per tu rb inepdf . i x [ i , ' t o t a l o p e x ' ] = OC. i x [ pdf . i x [ i , 'N ' ] , 'Cuml OC ' ] *opex per tu rb ineAEP = pdf . i x [ i , ' Power ' ] * 365 * 24 / 1E6pdf . i x [ i , 'AEP ' ] = AEPcapex = pdf . i x [ i , ' t o t a l c a p e x ' ]opex = pdf . i x [ i , ' t o t a l o p e x ' ]l i f e = f i n a n c i a l [ ' L i f e ' ]i r a t e = f i n a n c i a l [ ' i r a t e ' ]this LCOE = LCOE( capex , opex ,AEP, l i f e , i r a t e )pdf . i x [ i , 'LCOE ' ] = this LCOE# p r i n t LCOE(140000 ,0 ,35*365*24 ,10 ,0.15)i f i %1000 == 0:p r i n t ip r i n t ' [DONE] 'p r i n t ' P l o t t i n g . . . 'f i g u r e ( 1 )# p l o t ( pdf .AEP, pdf .LCOE)pdf .LCOE. p l o t ( logx=True )x l a b e l ( ' Number o f t u rb i nes ' )y l a b e l ( 'LCOE ( $ /MWh ' )f i g u r e ( 2 )pdf . p l o t ( 'AEP ' , 'LCOE ' , logx=True )x l a b e l ( ' Annual Energy Product ion (TWh/ y r ) ' )y l a b e l ( 'LCOE ( $ /MWh ' )o u t p u t l i s t . append ( pdf )p r i n t ' [DONE] 'pdf [ ' AEP per turb ine ' ] = pdf [ 'AEP ' ] / pdf [ 'N ' ]f . save f ig155f o r pdf i n o u t p u t l i s t :f1= f i g u r e ( 1 )pdf . p l o t ( 'N ' , 'LCOE ' , logx=True , s t y l e = ' k ' )x l a b e l ( ' Number o f i n s t a l l e d tu rb i nes ' )y l a b e l ( 'LCOE ( $ /MWh) ' )y l im ( [ 0 , 1 0 0 0 ] )x l im ( [1 ,13000 ] )# f1 . s e t s i z e i n c h e s ( 5 . 5 , 4 )f2= f i g u r e ( 2 )pdf . p l o t ( 'AEP ' , 'LCOE ' , logx=True , s t y l e = ' k ' )x l a b e l ( ' Annual Energy Product ion (MWh/ y r ) ' )y l a b e l ( 'LCOE ( $ /MWh) ' )x l im ( [ 1 E4,1E8 ] )y l im ( [ 0 , 1 0 0 0 ] )f3= f i g u r e ( 3 )pdf . p l o t ( 'N ' , 'AEP ' , s t y l e = ' k ' )y l a b e l ( 'AEP (MWhr/ y r ) ' )x l a b e l ( ' Number o f I n s t a l l e d Turbines ' )x l im ( [0 ,20000 ] )f4= f i g u r e ( 4 )pdf . p l o t ( 'N ' , ' AEP per turb ine ' , s t y l e = ' k ' , logy=True )x l a b e l ( ' Number o f I n s t a l l e d Turbines ' )y l a b e l ( 'AEP Per Turbine (MWhr/ y r ) ' )x l im ( [ 0 , 5 0 0 0 ] )f1 . save f ig ( 'LCOEvsN. eps ' )f2 . save f ig ( 'LCOEvsAEP. eps ' )f o r pdf i n o u t p u t l i s t :p r i n t pdf .LCOE. min ( )OC. p l o t ( 'N ' , 'OC ' , g r i d =False , s t y l e = ' k ' )x l im ( [ -50 ,1500 ] )y l im ( [ 0 , 1 ] )y l a b e l ( ' F rac t i on o f O r i g i n a l Cost ' )x l a b e l ( ' I n s t a l l e d Uni ts ' )save f ig ( ' l ea rn i ng . pdf ' )156Appendix CMatlab Code for Control ofLaboratory ExperimentC.1 ADV ProcessSerial.mf u n c t i o n arrayOut = ADV ProcessSer ia lSt r ing ( s t r I n )%For S t r i n g Inpu t%s t r I ntemparray = s t r s p l i t ( ' ' , s t r I n ) ;s i ze ( temparray )i f ( s i ze ( temparray == 10) )arrayOut = [ ] ;f o r i =1: leng th ( temparray )tempItem = st r2doub le ( temparray ( i ) ) ;i f ( isnan ( tempItem ) ~= 1)arrayOut =[ arrayOut , tempItem ] ;endendelsef p r i n t f ( ' Bad s ize ! ' )arrayOut = [ ]endC.2 ADV ProcessSerialString.mf u n c t i o n arrayOut = ADV ProcessSer ia lSt r ing ( s t r I n )157%For S t r i n g Inpu t%s t r I n ;temparray = s t r s p l i t ( ' ' , s t r I n ) ;arrayOut = [ ] ;f o r i =1: leng th ( temparray )tempItem = st r2doub le ( temparray ( i ) ) ;i f ( isnan ( tempItem ) ~= 1)arrayOut =[ arrayOut , tempItem ] ;endendi f ( s i ze ( arrayOut , 2 ) ~= 10)f p r i n t f ( ' bad s ize\n ' )arrayOut = [ ]endC.3 ADVclean.mf u n c t i o n out=ADVclean ( i n )th resho ld = 80; % %c o r r e l a t i o nCuto f fFreq = 1 / 3 ; %HzSampleFreq = 25;i f ( isnan ( i n ) )p r i n t f ( 'ERROR: ADV Data ar ray has a NaN. ' )r e t u r nendout = ADVclearLowCorrelat ion ( in , th resho ld ) ;i f ( l eng th ( out )>0)f o r i i =2:4out ( : , i i ) = f f t smoo th ( out ( : , i i ) , SampleFreq ,0 , Cuto f fFreq ) ;endendC.4 ADVclearLowCorrelation.mf u n c t i o n dataOut = ADVclearLowCorrelat ion ( dataIn , th resho ld )dataOut = [ ] ;f o r i =1: leng th ( data In )i f ( da ta In ( i , 8 )>t h resho ld )dataOut = [ dataOut ; data In ( i , : ) ] ;endendend158C.5 ADVgetdata.m%Get S e r i a l data from ADVADVdataBuffer = [ ] ;wh i le ( get ( ADVser ia lob ject , ' By tesAva i lab le ' )>100)ADVdataBuffer =[ ADVdataBuffer ; ADV ProcessSer ia lSt r ing ( f scan f ( ADVser ia lob jec t ) )] ;endC.6 ADVinit.mf p r i n t f ( ' I n i t i a l i z i n g ADV S e r i a l Communication . . . ' )ADVser ia lob jec t = s e r i a l ( 'COM1 ' , ' BaudRate ' , 19200) ;se t ( ADVser ia lob ject , ' I npu tBu f fe rS i ze ' ,4096)fopen ( ADVser ia lob jec t )C.7 ADVzero.m%ADVZerot i cADVZeroData = [ ] ;wh i le ( toc<120)ADVgetdata ;ADVZeroData = [ ADVZeroData ; ADVdataBuffer ] ;pause ( 0 . 0 5 ) ;i f (mod( toc , 5 )<0.1)f p r i n t f ( i n t 2 s t r ( toc ) )f p r i n t f ( '\n ' )endendADVZero = mean( ADVZeroData ( : , 2 ) )C.8 calbrateUSsensors.m%C a l i b r a t i o n Routine f o r US Sensorsmessages = { ' Sensor 1 , t a l l b lock ' , ' Sensor 1 , sho r t b lock ' , ' Sensor 2 , t a l l b lock ', ' Sensor 2 , shor t b lock ' , ' Sensor 3 , t a l l b lock ' , ' Sensor 3 , sho r t b lock ' , 'Sensor 4 , t a l l b lock ' , ' Sensor 4 , sho r t b lock ' } ;tempUScalib = UScal ib ;f o r i =1:8d isp ( messages ( i ) )159user in = i npu t ( ' Press c to c a l i b r a t e , s to sk ip : ' , ' s ' ) ;i f ( user in == ' s ' );e l s e i f ( user in == ' c ' )dataokay = 0;wh i le ( ~dataokay )s t a r t ( daqobject ) ;pause (10) ;stop ( daqobject ) ;[ tempdata , temptime ] = getdata ( daqobject , get ( daqobject , 'SamplesAvai lable ' ) ) ;p l o t ( tempdata ) ;user in2 = inpu t ( ' Data Okay? ( y / n ) : ' , ' s ' ) ;i f ( user in2 == ' y ' )f p r i n t f ( ' Data accepted . Moving on .\n ' ) ;dataokay = 1;e lsef p r i n t f ( ' Redoing measurement .\n ' ) ;endc lose a l l ;endtempUScalib ( i , 2 ) = mean( tempdata ( : , c e i l ( i / 2 ) ) ) ;endendC.9 calcPower.m%Calcu la te Power Outputt i cdel taH = u l t raC lean (WLs( : , 3 ) ) - us3zero - u l t raC lean ( (WLs( : , 2 ) ) - us2zero ) ;Q = data ( : , 2 ) -ADVZero ;%Assume deltaH has more po in t s than Q.Power = [ ] ;LQ = leng th (Q)LH = leng th ( del taH )tocExpQ = zeros (LH, 1 ) ;f o r i =1:LH%iQ = 1+ f l o o r (LQ* i / LH) )ExpQ( i ) = Q(1+ round ( ( LQ- 1 ) * i / LH) ) ;endPower = del taH . * ExpQ ;toc160C.10 CalcPowerErrorBars.mf u n c t i o n [ e r ro rbars , data , proper ] = CalcPowerErrorBars ( depthser ies , check )usSensorBound = 0.0015;number samples = 1E5 ;alpha = 2;i =1;us10 = P ( [ depthser ies ( i , 1 ) , depthser ies ( i , 1 ) -usSensorBound , depthser ies ( i , 1 ) +usSensorBound ] , alpha , number samples ) ;us20 = P ( [ depthser ies ( i , 2 ) , depthser ies ( i , 2 ) -usSensorBound , depthser ies ( i , 2 ) +usSensorBound ] , alpha , number samples ) ;us30 = P ( [ depthser ies ( i , 3 ) , depthser ies ( i , 3 ) -2* usSensorBound , depthser ies ( i , 3 ) +2*usSensorBound ] , alpha , number samples ) ;us40 = P ( [ depthser ies ( i , 4 ) , depthser ies ( i , 4 ) -usSensorBound , depthser ies ( i , 4 ) +usSensorBound ] , alpha , number samples ) ;% us10 = depthser ies ( i , 1 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us20 = depthser ies ( i , 2 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us30 = depthser ies ( i , 3 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us40 = depthser ies ( i , 4 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;Q0 = P ( [ check ( i , 1 ) , check ( i , 7 ) , check ( i , 8 ) ] , alpha , number samples ) ;%Q0 = check ( i , 7 ) + ( check ( i , 8 ) - check ( i , 7 ) ) . * rand ( number samples , 1 ) ;dt0 = us30 - us20 ;k = dt0 . / Q0. ˆ 2 ;ks i ng le = ( depthser ies ( i , 3 ) - depthser ies ( i , 2 ) ) / check (1 ,1 ) . ˆ 2proper= [ ] ;e r ro rba rs = [ ] ;f o r i =1: leng th ( depthser ies )% us1 = depthser ies ( i , 1 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us2 = depthser ies ( i , 2 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us3 = depthser ies ( i , 3 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;% us4 = depthser ies ( i , 4 ) + (2 * usSensorBound * rand ( number samples , 1 ) ) -usSensorBound ;us1 = P ( [ depthser ies ( i , 1 ) , depthser ies ( i , 1 ) *0 .99 - usSensorBound , depthser ies ( i, 1 ) *1.01+ usSensorBound ] , alpha , number samples ) ;us2 = P ( [ depthser ies ( i , 2 ) , depthser ies ( i , 2 ) *0 .99 -2* usSensorBound , depthser ies ( i, 2 ) *1.01+2* usSensorBound ] , alpha , number samples ) ;us3 = P ( [ depthser ies ( i , 3 ) , depthser ies ( i , 3 ) *0 .99 - usSensorBound , depthser ies ( i, 3 ) *1.01+ usSensorBound ] , alpha , number samples ) ;us4 = P ( [ depthser ies ( i , 4 ) , depthser ies ( i , 4 ) *0 .99 - usSensorBound , depthser ies ( i, 4 ) *1.01+ usSensorBound ] , alpha , number samples ) ;161% Q = check ( i , 7 ) + ( check ( i , 8 ) - check ( i , 7 ) ) . * rand ( number samples , 1 ) ;Q = P ( [ check ( i , 1 ) , check ( i , 7 ) , check ( i , 8 ) ] , alpha , number samples ) ;d t = us3 - us2 - k . *Q. ˆ 2 ;dh = us4 - us1 ;eta = dt . * Q . / (Q0 . * dh ) ;d t s i n g l e = ( depthser ies ( i , 3 ) - depthser ies ( i , 2 ) ) - ks i ng le * check ( i , 1 ) . ˆ 2 ;proper = [ proper ; d t s i n g l e * check ( i , 1 ) / ( ( depthser ies ( i , 4 ) - depthser ies ( i, 1 ) ) * check (1 ,1 ) ) ] ;data = [ us1 , us2 , us3 , us4 ,Q] ;e r ro rba rs = [ e r ro rba rs ; mean( eta ) - q u a n t i l e ( eta , 0 . 7 5 ) ] ;end%p l o t ( proper ) ;C.11 CalcShutterK.m%This i s l i k e experiment 5 but c o n t r o l l s f o r water l e v e l on the r i g h t s ide .%Update from Ex6 i n t h a t the r i g h t wei r i s t imed r a t h e r than using the%encoder .shu t te rSe t = [ 0 ]movement = 1000;%WLSet =[ 0 . 6 0 , 0 . 5 6 ; 0 . 5 9 , 0 . 5 6 ; 0 . 5 9 , 0 . 5 8 ; 0 . 5 9 , 0 . 5 7 ; 0 . 5 9 , 0 . 5 7 ; 0 . 5 8 , 0 . 5 7 ; 0 . 5 8 , 0 . 5 6 ; ] ; %mi f ( ~ e x i s t ( ' weirSP ' ) )weirSP = 0;endi f ( ~ e x i s t ( ' n ' ) )n = 0 ;endf o r WLi = [ 1 : 1 5 ]weir1move ( movement ) ;weir2move ( - movement ) ;t i cwh i le ( toc < 30)pause ( 0 . 5 ) ;i f (mod( toc , 5 )<0.49)toc162endendusdata = [ ] ;startDAQ ( ) ;t i cwh i le ( toc < 20)pause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;i f (mod( toc , 5 )<0.49)tocendendp r o f i l e Q% startDAQ ( ) ;n=n+1;ExData ( n ) . Data = usdataExData ( n ) . v e l o c i t i e s = Qdata ;ExData ( n ) . s h u t t e r = 0 ;ExData ( n ) . date = c lock ;c l k = c lockf i lename = [ date , ' ' , num2str ( c l k ( 4 ) ) , ' - ' , num2str ( c l k ( 5 ) ) ]save ( f i lename , ' ExData ' )endC.12 CleanUSspikes.mf u n c t i o n usdata = CleanUSspikes ( usdata )% subp lo t (2 ,1 ,1 ) ;% p l o t ( usdata ( : , 2 : 5 ) )f o r i i = 2:5%j k = f i n d ( usdata ( : , i )>4.9) ;j k = f i n d ( usdata ( : , i i )>(min ( usdata ( : , i i ) + .8 ) ) ) ;i f ( ~ isempty ( j k ) )k i l l I n d i c e s = [ ] ;f o r i =1: leng th ( j k )k i l l I n d i c e s = [ k i l l I n d i c e s , ( j k ( i ) -20) : ( j k ( i ) +20) ] ;endk i l l I n d i c e s = k i l l I n d i c e s ( f i n d ( k i l l I n d i c e s >=1) ) ;k i l l I n d i c e s = k i l l I n d i c e s ( f i n d ( k i l l I n d i c e s <=leng th ( usdata ) ) ) ;usdata ( k i l l I n d i c e s , i i ) =nan ;end163end% subp lo t (2 ,1 ,2 ) ;% p l o t ( usdata ( : , 2 : 5 ) ) ;% drawnow ;% pause ( 2 ) ;C.13 countdown.mf u n c t i o n countdown ( t ime )t i c ;wh i le ( toc<t ime )i f (mod( toc , 5 )<0.001)f p r i n t f ( i n t 2 s t r ( t ime - toc ) )f p r i n t f ( ' . ' )endendf p r i n t f ( '\n ' ) ;C.14 daqinit.mdaqobject = ana log input ( 'mcc ' ,0 ) ;daq sample rate = 30;%Analog Inpu t Parametersaddchannel ( daqobject , 0 : 3 ) ;se t ( daqobject , ' SampleRate ' , daq sample rate ) ;se t ( daqobject , ' SamplesPerTrigger ' ,10000) ;se t ( daqobject , ' Tr iggerRepeat ' , i n f ) ;C.15 define global.mg loba l s h u t t e r S e r i a l O b j e c t ;g loba l we i r1Ser ia lOb jec t ;g loba l we i r2Ser ia lOb jec t ;g loba l weir1zero ;g loba l daqobject ;g loba l ADVser ia lob jec t ;g loba l t r a v S e r i a l O b j e c t ;g loba l UScal ibrho = 998;g = 9 .81 ;164C.16 DetermineFlow.mf u n c t i o n [Q, uMatr ix , Qmin , Qmax] = DetermineFlow ( ve lS t ruc t , t o t a l d e p t h )%Gets passed ExData . v e l o c i t i e s s t r u c t u r edepth from bot tom = 0.028;steps per m = 35000 / .278 ;ve l sca le = 1/10000; %Converts ADV v e l o c i t y to m/ schannel breadth = 28.25*0.0254;ys = [ ] ;xs = [ ] ;f o r profN = 1: leng th ( v e l S t r u c t ) %Loops through each measurement l o c a t i o nys = [ ys ; v e l S t r u c t ( profN ) . y ] ;xs = [ xs ; v e l S t r u c t ( profN ) . x ] ;end%Co l l ec t s a l l the unique l o c a t i o n s o f measurements f o r x and y ( there are%doubles , otherwise )yset = unique ( ys ) ;xset = unique ( xs ) ;uMatr ix = zeros ( leng th ( xset ) , l eng th ( yset ) ) ;f o r profN = 1: leng th ( v e l S t r u c t ) %Loops through each measurement l o c a t i o nxindex = f i n d ( xset == v e l S t r u c t ( profN ) . x ) ;y index = f i n d ( yset == v e l S t r u c t ( profN ) . y ) ;VelData = ADVclean ( v e l S t r u c t ( profN ) . u ) ;%uMatr ix ( xindex , y index ) = mean( norm ( [ VelData ( : , 2 ) , VelData ( : , 3 ) , VelData ( : , 4 ) ] )) ;uMatr ix ( xindex , y index ) = mean ( ( VelData ( : , 2 ) . ˆ2+ VelData ( : , 3 ) . ˆ2+ VelData ( : , 4 ). ˆ 2 ) . ˆ 0 . 5 ) ;enduMatr ix ;%Convert the uMatr ix i n t o a f low ra te%Calcu la te the s p e c i f i c f l ow f o r each v e r t i c a l s l i c eqset = [ ] ;f o r x i = 1 : leng th ( xset )q = 0;f o r y i = 1 : leng th ( yset )%Note : area i s s p e c i f i c - area ( i e length , not area ) . Deal t w i th%l a t e r .i f ( y i == 1)%Bottom measurementArea = ( yset ( y i ) / steps per m ) +depth from bot tom ;ve l = ve l sca le * uMatr ix ( x i , y i ) / 2 ; %Assumes no s l i p on the bottomq = q + Area * ve l ;e l s e i f ( y i == leng th ( yset ) )165%Top measurementArea = t o t a l d e p t h - ( ( yset ( y i ) / steps per m ) +depth from bot tom ) ;ve l = ve l sca le * uMatr ix ( x i , y i ) ; %Assumes no s l i p on the bottomq = q + Area * ve l ;e lse%a l l the middle onesArea = ( yset ( y i ) - yset ( y i - 1 ) ) / steps per m ;ve l = ve l sca le * ( uMatr ix ( x i , y i ) +uMatr ix ( x i , y i - 1 ) ) / 2 ; %Assumes no s l i pon the bottomq = q + Area * ve l ;endendqset = [ qset ; q ] ;endi f ( l eng th ( xset ) == 3)%Q = channel breadth * ( qset ( 1 ) +2* qset ( 2 ) +qset ( 3 ) ) / 4 ;Q = channel breadth *mean( qset ) ;Qmin = channel breadth * min ( qset ) * . 9 ;Qmax = channel breadth *max( qset ) * 1 . 1 ;% Q = channel breadth * qset ( 2 ) ;e l s e i f ( l eng th ( xset ) == 4)Q = channel breadth *mean( qset ) ;%Q = channel breadth * ( qset ( 1 ) +2* qset ( 2 ) +2* qset ( 3 ) +qset ( 4 ) ) / 6 ;Qmin = channel breadth * min ( qset ) * . 9 ;Qmax = channel breadth *max( qset ) * 1 . 1 ;e l s e i f ( l eng th ( xset ) == 1)Q = channel breadth * qset ( 1 ) ;Qmin = channel breadth * min ( qset ) * . 9 ;Qmax = channel breadth *max( qset ) * 1 . 1 ;e lsexsete r r o r ( ' Unknown xset s ize ! ' ) ;endC.17 DispDepths.m%Reads and output depths and dHdepths = [ ] ;wh i le ( 1 )usdata = [ ] ;startDAQ ( ) ;t i cwh i le ( toc < 2) %C o l l e c t data f o r 2pause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;% i f (mod( toc , 5 )<0.49)% toc% endend166%El im ina te a l l the da tapo in ts where the sensor looses ' lock ' .usdata = CleanUSspikes ( usdata ) ;% p l o t ( usdata ( : , 2 ) ) ;RightWL = usdepth (mynanmean( usdata ( : , 2 ) ) ,1 ) %Find mean water l e v e lLeftWL = usdepth (mynanmean( usdata ( : , 5 ) ) ,4 ) %Find mean water l e v e lMidLeftWL = usdepth (mynanmean( usdata ( : , 4 ) ) ,3 ) %Find mean water l e v e lMidRightWL = usdepth (mynanmean( usdata ( : , 3 ) ) ,2 ) %Find mean water l e v e ldepths = [ depths ; RightWL , MidRightWL , MidLeftWL , LeftWL ] ;hold on ;p l o t ( depths ( : , 1 ) , ' b ' ) ;p l o t ( depths ( : , 2 ) , ' g ' ) ;p l o t ( depths ( : , 3 ) , ' r ' ) ;p l o t ( depths ( : , 4 ) , ' c ' ) ;drawnow ;endC.18 Ex7.m%This i s l i k e experiment 5 but c o n t r o l l s f o r water l e v e l on the r i g h t s ide .%Update from Ex6 i n t h a t the r i g h t wei r i s t imed r a t h e r than using the%encoder .shu t te rSe t = [28E5:4E5:4E6 ]WLSet = [0 .545 , . 524 ]%WLSet =[ 0 . 6 0 , 0 . 5 6 ; 0 . 5 9 , 0 . 5 6 ; 0 . 5 9 , 0 . 5 8 ; 0 . 5 9 , 0 . 5 7 ; 0 . 5 9 , 0 . 5 7 ; 0 . 5 8 , 0 . 5 7 ; 0 . 5 8 , 0 . 5 6 ; ] ; %mi f ( ~ e x i s t ( ' weirSP ' ) )weirSP = 0;endi f ( ~ e x i s t ( ' n ' ) )n = 0 ;endf o r WLi = 1 : leng th (WLSet )LeftWL SP = WLSet ( WLi , 1 ) ;RightWL SP = WLSet ( WLi , 2 ) ;weir1move ( 0 ) ;weirSP = 0;f o r s h u t t e r = shu t te rSe tSeeking = [ ] ;shuttermove ( s h u t t e r )t i cwr ( ' Wai t ing . . . ' ) ;wh i le ( toc < 20)i f (mod( toc , 5 )<0.49)tocend167pause ( . 5 ) ;end% PIDWeirs ;WLthreshold = 0.00025;RightWLerror = WLthreshold *10 ; %Force loop .Lef tWLerror = WLthreshold *10 ; %Force loop .%I t e r a t e to set l e f t water l e v e lwh i le ( ( abs ( Lef tWLerror ) > WLthreshold ) | ( abs ( RightWLerror ) > WLthreshold ) )%C o l l e c t US data to f i n d cu r ren t water l e v e l a t Right wei rusdata = [ ] ;startDAQ ( ) ;t i cwh i le ( toc < 10) %C o l l e c t data f o r 10spause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;i f (mod( toc , 5 )<0.49)tocendend%El im ina te a l l the da tapo in ts where the sensor looses ' lock ' .usdata = CleanUSspikes ( usdata ) ;% p l o t ( usdata ( : , 2 ) ) ;RightWL = usdepth (mynanmean( usdata ( : , 2 ) ) ,1 ) %Find mean water l e v e lRightWLerror = RightWL - RightWL SP %Find e r r o r- RightWLerror *3*62.5*1000%Actuate the wei ri f ( abs ( RightWLerror ) > WLthreshold )weir1move ( round ( - RightWLerror *2*62.5*1000) ) ;%Wait f o r 20 seconds .t i cwr ( ' Wai t ing f o r adjustment ' ) ;Lef tWLerror = WLthreshold *15 ;wh i le ( toc < 20)pause ( 0 . 5 ) ;i f (mod( toc , 5 )<0.49)tocendend %End of wa i tend%%%Now actuate the l e f t wei r%C o l l e c t US data to f i n d cu r ren t water l e v e l a t l e f t wei r168usdata = [ ] ;startDAQ ( ) ;t i cwh i le ( toc < 10) %C o l l e c t data f o r 10spause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;i f (mod( toc , 5 )<0.49)tocendend%El im ina te a l l the da tapo in ts where the sensor looses ' lock ' .usdata = CleanUSspikes ( usdata ) ;% p l o t ( usdata ( : , 2 ) ) ;LeftWL = usdepth (mynanmean( usdata ( : , 5 ) ) ,4 ) %Find mean water l e v e lLef tWLerror = LeftWL - LeftWL SP %Find e r r o r- Lef tWLerror *3*62.5*1000%Actuate the wei ri f ( abs ( Lef tWLerror ) > WLthreshold )weir2move ( round ( - Lef tWLerror *2*62.5*1000) ) ;%Wait f o r 30 seconds .t i cwr ( ' Wai t ing f o r L e f t adjustment ' ) ;RightWLerror = WLthreshold *10 ;wh i le ( toc < 30)pause ( 0 . 1 ) ;i f (mod( toc , 5 )<0.49)tocpause ( 0 . 5 )endend %End of wa i tendSeeking = [ Seeking ; RightWL , LeftWL ]hold on ;p l o t ( Seeking ( : , 1 ) , ' r ' ) ;p l o t ( Seeking ( : , 2 ) , ' b ' ) ;drawnow ;end %end of wei r loopt i cusdata = [ ] ;startDAQ ( ) ;t i cwh i le ( toc < 20)pause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;169i f (mod( toc , 5 )<0.49)tocendendp r o f i l e Q% startDAQ ( ) ;t i cwh i le ( toc < 20)pause ( 0 . 5 ) ;usdata = [ usdata ; newDAQdata ( ) ] ;i f (mod( toc , 5 )<0.49)tocendendn=n+1;ExData ( n ) . Data = usdataExData ( n ) . v e l o c i t i e s = Qdata ;ExData ( n ) . weir1 = weirSP ;ExData ( n ) . s h u t t e r = s h u t t e r ;ExData ( n ) . date = c lock ;c l k = c lockExData ( n ) . RightWL = RightWL SP ;ExData ( n ) . LeftWL = LeftWL SP ;ExData ( n ) . Seeking = Seeking ;ExData ( n ) . USCal ib ra t ion = UScal ib ;f i lename = [ date , ' ' , num2str ( c l k ( 4 ) ) , ' - ' , num2str ( c l k ( 5 ) ) ]save ( [ ' data / ' , f i lename ] , ' ExData ' )endendC.19 fftsmooth.mf u n c t i o n y= f f t smoo th ( y , fs , lowband , highband )%FFT Smooth% nChan i s 2 f o r a s tereo f i l e : Independent l e f t / r i g h t channeln = leng th ( y ) ;% c a l c u l a t e the frequency corresponding to each FFT bin .% t h i s inc ludes negat ive f requenc ies !f reqbase= f s * (mod ( ( ( 0 : n - 1 ) + f l o o r ( n / 2 ) ) , n ) - f l o o r ( n / 2 ) ) / n ;% b u i l d a mask [ . . . 0 0 0 1 1 . . . 1 1 0 0 0 . . . ]% wi th one en t ry per frequencymask=( abs ( freqbase ) > lowband ) . * ( abs ( freqbase ) < highband ) ;170mask=mask ' ; % conver t to column vec to r ( the audio s i g n a l i s a lso a column vec to r )% e x t r a c t the channelyc=y ;% to frequency domainyc= f f t ( yc ) ;% f i g u r e ( )% p l o t ( freqbase , abs ( yc ) )% apply f i l t e ryc = yc . * mask ;% back to t ime domainyc= i f f t ( yc ) ;% ove rwr i t e o r i g i n a l data% r e a l ( ) i s needed to remove imaginary pa r t caused by f i n i t e numer icalp r e c i s i o ny ( : ) = r e a l ( yc ) + mean( y ) ;% In extreme cases , r i n g i n g can increase the h ighes t peaks .% Therefore , normal ize the peaks to 1% peak=max(max( abs ( y ) ) ) ;% y=y * 1 / peak ;% w r i t e output f i l eC.20 int2Arduino.mf u n c t i o n s t r i ngOu t = in t2Ardu ino ( pos )d i g i t L e n g t h = 8;i f abs ( pos )<10ˆ d i g i t L e n g t hs t r i ngOu t = ' saaa ' ;i f ( pos < 0)middle = ' ccc ' ;pos = -pos ;e lsemiddle = ' bbb ' ;endnumString = zeros (1 , d i g i t L e n g t h ) ;f o r i = (1 : d i g i t L e n g t h )pow = d i g i t L e n g t h - i +1;numString (pow) = mod( pos ,10 ) + ' 0 ' ;pos = f i x ( pos /10 ) ;ends t r i ngOu t = [ ' saaa ' , numString , middle , numString , ' xxx ' ] ;171elses t r i ngOu t = '\n ' ;endendC.21 int2string.mf u n c t i o n s t r i n g = i n t 2 s t r i n g ( n )n = f i x ( n ) ;i f ( n<0)s t r i n g = ' - ' ;n = abs ( n ) ;e lses t r i n g = ' ' ;endi n t A r r a y = [ ] ;wh i le ( n>0)i n t A r r a y = [ i n tA r ray ,mod( n ,10 ) ] ;n = f i x ( n /10 ) ;endf o r i = 1 : leng th ( i n t A r r a y )s t r i n g = [ s t r i n g , i n t A r r a y ( leng th ( i n t A r r a y ) - i +1)+ ' 0 ' ] ;endendC.22 mergeDAQ.mf u n c t i o n expADV=mergeDAQ( data ,ADV, ADVRate )t i c%Gets mat r i x o f data from the analog to d i g i t a l board and i n s e r t s ADV%measurements w i th i t .%Data ar ray column 1 i s the t ime i n seconds .%ADV column 1 i s the sample number , not re fe r renced to zero .%The beginnign o f both m a t r i c i e s occurs a t the same po in t .%DAQ Data should be at a h igher sampling ra te than the ADV.ADVtime = (ADV( : , 1 ) -ADV(1 ,1 ) ) / ADVRate ;expADV = [ ] ;f o r i =2:4expADVx = [ ] ;%Purposebu i l t a lgo r i t hmf o r j j =2: leng th (ADV)ind i ces = f i n d ( ( data ( : , 1 )>=ADVtime ( j j - 1 ) ) &( data ( : , 1 )<ADVtime ( j j ) ) ) ;i f ( l eng th ( i nd i ces )>0)%mult = ( 0 : 1 / leng th ( i nd i ces ) : ( 1 - 1 / leng th ( i nd i ces ) ) ) ' ;172expADVx = [ expADVx ; ( ( ADV( j j , i ) -ADV( j j -1 , i ) ) / ( ADVtime ( j j ) -ADVtime ( j j- 1 ) ) ) . * ( data ( ind ices , 1 ) -ADVtime ( j j - 1 ) ) + ADV( j j -1 , i ) ] ;%expADVx = [ expADVx ; ADV( j j -1 , i ) + mult . * (ADV( j j , i ) - ADV( j j -1 , i ) )] ;endendexpADV = [ expADV , expADVx ] ;% Code using b u i l t i n i n t e r p o l a t i o n f u n c t i o n . Slow .% f o r i i =1: leng th ( data )% expADV( i i ) = [ i n t e r p 1 ( ADVtime ,ADV( : , i ) , data ( i i , 1 ) ) ] ;% endtocendC.23 mynanmean.mf u n c t i o n x = mynanmean( i n )inds = f i n d ( ~ isnan ( i n ) ) ;x = mean( i n ( inds ) ) ;C.24 newADVdata.mf u n c t i o n ADVdataBuffer=newADVdata ( )g loba l ADVser ia lob jec t ;ADVdataBuffer = [ ] ;wh i le ( get ( ADVser ia lob ject , ' By tesAva i lab le ' )>100)ADVdataBuffer =[ ADVdataBuffer ; ADV ProcessSer ia lSt r ing ( f scan f ( ADVser ia lob jec t ) )] ;endC.25 newDAQdata.mf u n c t i o n data = newDAQdata ( )g loba l daqobject ;data = [ ] ;i f ( get ( daqobject , ' SamplesAvai lable ' )>0)%SampsAvail = get ( daqobject , ' SamplesAvai lable ' ) ;[ tempdata , temptime ] = getdata ( daqobject , get ( daqobject , ' SamplesAvai lable ' ) ) ;data = [ temptime , tempdata ] ;e lser e t u r nend173C.26 ProcessEx5.mFlumeWidth = 0 . 5 ; %mExData2 = s t r u c t ( ' meanDepth ' ,{} , 'Q ' ,{} , ' Power ' ,{} , ' Shu t te r ' ,{} , 'dH ' ,{} , ' dT ' ,{} , 'depth ' ,{} )check = [ ]i f ( ~ e x i s t ( ' Nset ' ) )Nset = 1 : leng th ( ExData ) ;enddepthser ies = [ ] ;f o r exN = NsetexN ;ADV1 = [ ] ;f o r profN = 1: leng th ( ExData (exN) . v e l o c i t i e s )i f ( ~ isempty ( ExData (exN) . v e l o c i t i e s ( profN ) . u ) )ADV1 = [ADV1; ADVclean ( ExData (exN) . v e l o c i t i e s ( profN ) . u ) ] ;endendV = mean(ADV1( : , 2 ) ) /10000; %m/ sUS1 = CleanUSspikes ( ExData (exN) . Data ) ;depth = [ ] ;meanDepth = [ ] ;f o r usN = 2:5depth = [ depth , USdepth (US1 ( : , usN) , ( usN - 1 ) ) ] ;meanDepth = [ meanDepth ; nanmean ( depth ( : , usN - 1 ) ) ] ;endExData2 (exN) . depth = depth ;depthser ies = [ depthser ies ; meanDepth ' ] ;ExData2 (exN) . meanDepth = meanDepth ;ExData2 (exN) . Shut te r = ExData (exN) . s h u t t e r ;[Q, uMatr ix , Qmin , Qmax] = DetermineFlow ( ExData (exN) . v e l o c i t i e s , ExData2 (exN) .meanDepth ( 3 ) ) ;ExData2 (exN) .Q = Q;ExData2 (exN) . Qmin = Qmin ;ExData2 (exN) .Qmax = Qmax;ExData2 (exN) .Q2 = FlumeWidth * ExData2 (exN) . meanDepth ( 3 ) * V ;ExData2 (exN) . dT = ExData2 (exN) . meanDepth ( 3 ) -ExData2 (exN) . meanDepth ( 2 ) ;ExData2 (exN) .dH = ExData2 (exN) . meanDepth ( 4 ) -ExData2 (exN) . meanDepth ( 1 ) ;ExData2 (exN) . uMat r ix = uMatr ix ;check = [ check ; ExData2 (exN) .Q, ExData2 (exN) . dT , ExData2 (exN) . dH,V, ExData2 (exN). Shut ter , ExData2 (exN) . dT* ExData2 (exN) .Q, Qmin ,Qmax ] ;end174C.27 ProcessMultipleExperiments.mg loba l UScal ib ;UScal ib = load ( ' USca l ib ra t ion -2011 -10 -28. t x t ' ) ;rho = 998;g=9.81;epr = s t r u c t ( ' d i r e c t o r y ' ,{} , ' f i lename ' ,{} , ' nmin ' ,{} , ' nmax ' ,{} , ' depths ' ,{} , 'da tapo in ts ' ,{} , ' U S c a l i b r a t i o n f i l e ' ,{} , ' DTca l i b ra t i on ' ,{} ) ;f i l e n = 0;f i l e n = f i l e n +1;epr ( f i l e n ) . d i r e c t o r y = ' data ' ;epr ( f i l e n ) . f i lename = ' 26 -Oct -2011 23 -36 . mat ' ;f i l e n = f i l e n +1;epr ( f i l e n ) . d i r e c t o r y = ' data ' ;epr ( f i l e n ) . f i lename = ' 29 -Oct -2011 18 -35 . mat ' ;epr ( f i l e n ) . nmin = 19;f i l e n = f i l e n +1;epr ( f i l e n ) . d i r e c t o r y = ' data ' ;epr ( f i l e n ) . f i lename = ' 04 -Nov-2011 0 -31 . mat ' ;% f i l e n = f i l e n +1;% epr ( f i l e n ) . d i r e c t o r y = ' data ' ;% epr ( f i l e n ) . f i lename = '26 - Oct -2011 16 -54 . mat ' ;% epr ( f i l e n ) . nmin = 14;%% f i l e n = f i l e n +1;% epr ( f i l e n ) . d i r e c t o r y = ' data ' ;% epr ( f i l e n ) . f i lename = '21 - Oct -2011 0 -50 . mat ' ;% f i l e n = f i l e n +1;% epr ( f i l e n ) . d i r e c t o r y = ' data ' ;% epr ( f i l e n ) . f i lename = '14 - Oct -2011 17 -43 . mat ' ;c lose a l lf o r fn =1: f i l e nload ( [ epr ( fn ) . d i r e c t o r y , ' / ' , epr ( fn ) . f i lename ] ) ;%Def ine the da tapo in ts to process - there may be extraneousi f ( isempty ( epr ( fn ) . nmin ) )nmin = 1;e lsenmin = epr ( fn ) . nmin ;end175i f ( e x i s t ( ' epr ( fn ) . USCal ib ra t ion ' ) )i f ( isempty ( epr ( fn ) . USCal ib ra t ion ) )UScal ib = load ( ' USca l ib ra t ion -2011 -10 -28. t x t ' ) ;e lseUScal ib = epr ( fn ) . USCal ib ra t ion ;endelseUScal ib = load ( ' USca l ib ra t ion -2011 -10 -28. t x t ' ) ;endi f ( 1 )UScal ib ( 1 : 2 , 4 ) = 11;UScal ib ( 3 : 4 , 4 ) = 11.6UScal ib ( 7 : 8 , 4 ) = 12 .5 ;endi f ( isempty ( epr ( fn ) . nmax) )nmax = leng th ( ExData ) ;e lsenmax = epr ( fn ) . nmax ;endNset = nmin : nmax ;%Process the data - t h i s creates the ' check ' mat r i xProcessEx5epr ( fn ) . da tapo in ts = check ;epr ( fn ) . depths = depthser ies ;%DT Cor rec t ion% 1 - quadra t i c based on f i r s t po i n t o f d t% 2 - p o l y f i t based from other experiment% 3 - d i r e c t i n t e r p o l a t i o n from preveous experimentmode = 1;swi tch modecase 1k = check (1 ,2 ) / check (1 ,1 ) ˆ2epr ( fn ) . QuadCorrectedDT = check ( : , 2 ) - k . * check ( : , 1 ) . ˆ 2 ;d t = epr ( fn ) . QuadCorrectedDT ;case 2dt = epr ( fn ) . PolyCorDT ;case 3epr ( fn ) . InterpCorDT = check ( : , 2 ) - i n t e r p 1 (nomDT( : , 1 ) ,nomDT( : , 2 ) ,check ( : , 1 ) ) ; %Corrected d t d i r e c t l y from datad t = epr ( fn ) . InterpCorDT ;otherwised isp ( ' E r ro r - need appropr ia te mode ' ) ;endDT In te rpo la ted = 0;i f ( DT In te rpo la ted )e lse176endepr ( fn ) . power = dt . * epr ( fn ) . da tapo in ts ( : , 1 ) * rho *g ;epr ( fn ) . k i = check (1 ,3 ) / check (1 ,1 ) ˆ 2 ;epr ( fn ) . d t = d t ;epr ( fn ) . k t = d t . / check ( : , 1 ) . ˆ 2 ;epr ( fn ) . Per rorbars = CalcPowerErrorBars ( epr ( fn ) . depths , epr ( fn ) . da tapo in ts ) ;epr ( fn ) . the ta = ( as in ( ( epr ( fn ) . da tapo in ts ( : , 5 ) -2 .2E6) . /3 .11127E6) ) *180/ p i +45;end%%plot num = 3close a l lf = f i g u r e ( ) ;hold on ;y l im ( [ 0 . 5 1 5 , 0 . 5 5 ] )p l o t ( epr ( plot num ) . theta , epr ( plot num ) . depths ( : , 1 ) , ' k ' ) ;p l o t ( epr ( plot num ) . theta , epr ( plot num ) . depths ( : , 2 ) , ' k - - ' ) ;p l o t ( epr ( plot num ) . theta , epr ( plot num ) . depths ( : , 3 ) , ' k - - ' ) ;p l o t ( epr ( plot num ) . theta , epr ( plot num ) . depths ( : , 4 ) , ' k ' ) ;se t ( gca , ' YTick ' , [ 0 . 5 2 : 0 . 0 1 : . 5 5 ] )se t ( gca , ' FontSize ' ,12)x l a b e l ( ' Shu t te r Angle ( degrees ) ' ) ;y l a b e l ( ' Water Surface E leva t ion (m) ' ) ;f i g u r e ( ) ;subp lo t (2 ,1 ,1 ) ;hold on ;[ haxes , h l ine1 , h l i ne2 ] = p l o t y y ( epr ( plot num ) . theta , epr ( plot num ) . dt , epr ( plot num ). theta , epr ( plot num ) . da tapo in ts ( : , 1 ) *1000) ;se t ( h l ine1 , ' L ineSty le ' , ' - ' )se t ( h l ine1 , ' Color ' , ' k ' )se t ( h l ine1 , ' Marker ' , ' * ' )se t ( h l ine2 , ' L ineSty le ' , ' - - ' )se t ( h l ine2 , ' Color ' , ' k ' )se t ( h l ine2 , ' Marker ' , ' o ' )se t ( gca , ' YColor ' , ' k ' )se t ( haxes ( 1 ) , ' y l im ' , [ - 0 . 0 0 4 , . 0 3 ] )se t ( haxes ( 1 ) , ' YTick ' , [ 0 : 0 . 0 0 5 : 0 . 0 3 ] )se t ( haxes ( 1 ) , ' YColor ' , ' k ' )se t ( haxes ( 1 ) , ' Box ' , ' o f f ' )se t ( haxes ( 2 ) , ' YTick ' , [ 0 : 2 5 : 1 5 0 ] )se t ( haxes ( 2 ) , ' YColor ' , ' k ' )axes ( haxes ( 1 ) )se t ( gca , ' FontSize ' ,12)y l a b e l ( ' Turbine Head h T (m) ' )177axes ( haxes ( 2 ) )se t ( gca , ' FontSize ' ,12)y l a b e l ( ' Flow Rate ( L / s ) ' )legend ( 'Q ' , ' h t ' )x l a b e l ( ' ( a ) Shut te r Angle ( degrees ) ' ) ;subp lo t (2 ,1 ,2 )p l o t ( epr ( plot num ) . theta , epr ( plot num ) . power , ' k - ˆ ' )se t ( gca , ' FontSize ' ,12)x l a b e l ( ' ( b ) Shut te r Angle ( degrees ) ' ) ;y l a b e l ( ' Ex t rac ted Power (W) ' ) ;se t ( gca , ' Box ' , ' o f f ' )f i g u r e ( )se t ( gca , ' FontSize ' ,11)hold on ;tempn= 1p l o t ( epr ( tempn ) . da tapo in ts ( : , 1 ) . / epr ( tempn ) . da tapo in ts (1 ,1 ) , epr ( tempn ) . power . / (epr ( tempn ) . da tapo in ts (1 ,1 ) . * epr ( tempn ) . da tapo in ts ( : , 3 ) * rho *g ) , ' k * ' )tempn= 2p l o t ( epr ( tempn ) . da tapo in ts ( : , 1 ) . / epr ( tempn ) . da tapo in ts (1 ,1 ) , epr ( tempn ) . power . / (epr ( tempn ) . da tapo in ts (1 ,1 ) . * epr ( tempn ) . da tapo in ts ( : , 3 ) * rho *g ) , ' ko ' )tempn= 3p l o t ( epr ( tempn ) . da tapo in ts ( : , 1 ) . / epr ( tempn ) . da tapo in ts (1 ,1 ) , epr ( tempn ) . power . / (epr ( tempn ) . da tapo in ts (1 ,1 ) . * epr ( tempn ) . da tapo in ts ( : , 3 ) * rho *g ) , ' k ˆ ' )se t ( gca , ' Xd i r ' , ' reverse ' )QonQ = 0 : 0 . 0 1 : 1 ;eta = QonQ-QonQ. ˆ 3 ;p l o t (QonQ, eta , ' k - - ' ) ;x l im ( [ . 2 5 , 1 . 0 5 ] )y l im ( [ - 0 . 0 2 , . 5 ] ) ;y l a b e l ( ' E t r a c t i o n Rat io \eta ' ) ;x l a b e l ( ' Flow ra te as a r a t i o o f n a t u r a l (Q/ Q 0 ) ' ) ;f i g u r e ( )%hard coded - need to f i xyset = [ 0;7500;15000;22500;30000;37500;45000;52500;60000]depth f rom bot tom = 0.028;steps per m = 35000 / .278 ;yset = yset / steps per m + depth from bot tom ;f i g u r e ( )hold a l lu P r o f i l e = [ ]f o r i j = 1 : leng th ( ExData2 )f o r j k = 1 : s ize ( ExData2 ( i j ) . uMatr ix , 1 )temp = ExData2 ( i j ) . uMatr ix ( j k , : ) ;u P r o f i l e = [ u P r o f i l e , temp ' / max( temp ) ] ;p l o t ( temp ' / max( temp ) , yset )endend178y l im ( [ 0 , 0 . 5 5 ] ) ;x l a b e l ( ' Normalized V e l o c i t y ( u / u {max}) ' )y l a b e l ( ' Depth (m) ' ) ;f i g u r e ( )se t ( gca , ' FontSize ' ,12)hold on ;fn= 1e r ro rba r ( epr ( fn ) . k t / epr ( fn ) . k i , epr ( fn ) . power . / ( epr ( fn ) . da tapo in ts ( : , 3 ) . * epr( fn ) . da tapo in ts (1 ,1 ) * rho *g ) , epr ( fn ) . Perrorbars , ' k * ' )fn= 2e r ro rba r ( epr ( fn ) . k t / epr ( fn ) . k i , epr ( fn ) . power . / ( epr ( fn ) . da tapo in ts ( : , 3 ) . * epr( fn ) . da tapo in ts (1 ,1 ) * rho *g ) , epr ( fn ) . Perrorbars , ' ko ' )fn= 3e r ro rba r ( epr ( fn ) . k t / epr ( fn ) . k i , epr ( fn ) . power . / ( epr ( fn ) . da tapo in ts ( : , 3 ) . * epr( fn ) . da tapo in ts (1 ,1 ) * rho *g ) , epr ( fn ) . Perrorbars , ' k ˆ ' )k t k i = 0 : 0 . 0 1 : 1 5 ;eta = k t k i . / ( k t k i + 1) . ˆ 1 . 5 ;p l o t ( k t k i , eta ) ;y l im ( [ 0 , 0 . 6 ] ) ;x l im ( [ - 0 . 0 1 , 1 1 ] ) ;y l a b e l ( ' E x t r a c t i o n Rat io \eta ' ) ;x l a b e l ( ' Re la t i ve Resistance k t / k i ' ) ;f i g u r e ( )hold a l l ;f o r fn =1: f i l e np l o t ( epr ( fn ) . k t / epr ( fn ) . k i , epr ( fn ) . da tapo in ts ( : , 1 ) . / epr ( fn ) . da tapo in ts (1 ,1 ) , '. ' )endQonQ0 = 1 . / ( k t k i +1) . ˆ 0 . 5 ;p l o t ( k t k i ,QonQ0) ;%% f i g u r e ( )% hold a l l ;% f o r fn =1: f i l e n% subp lo t ( f i l e n ,1 , fn )% p l o t ( epr ( fn ) . depths )% end%% f i g u r e ( )%% f o r fn =1: f i l e n% subp lo t ( f i l e n ,1 , fn ) ;% hold a l l ;% p l o t ( epr ( fn ) . da tapo in ts ( : , 2 ) ) ;% p l o t ( d t ) ;% end%coefs i s the f i t t i n g c o e f f i c i e n c e s from determin ing the gate h t to Q179%responseC.28 ProcessProfile.mf u n c t i o n [ u , v ,w,Q, x , y ]= ProcessPro f i l e ( ExID , ExTime )%Loads the f i l e s from a p r o f i l i n g run , c leans and processes .%Load the parametersf o l d e r = ' QPro f i l es / ' ;f i lename = [ fo l de r , 'E ' , num2str ( ExID ) , ' /E ' , num2str ( ExID ) , ' -T ' , num2str ( ExTime ) , ' .mat ' ] ;load ( f i lename ) ;%Calcu la te the p o s i t i o n sx i = 1 : num width samples ;y i = 1 : num height samples ;x = ( x i - 1 ) * max width / ( num width samples - 1 ) ;y = ( y i - 1 ) * max height / ( num height samples - 1 ) ;u=zeros ( leng th ( x ) , leng th ( y ) ) ;v=zeros ( leng th ( x ) , leng th ( y ) ) ;w=zeros ( leng th ( x ) , leng th ( y ) ) ;f o r i i = 1 : leng th ( x )f o r j j = 1 : leng th ( y )x i = x ( i i ) ;y i = y ( j j ) ;f i lename = [ fo l de r , 'E ' , num2str ( ExID ) , ' /E ' , num2str ( ExID ) , ' -T ' , num2str (ExTime ) , ' -X ' , num2str ( x i ) , ' -Y ' , num2str ( y i ) , ' . t x t ' ]data = load ( f i lename ) ;data = ADVclean ( data ) ;u ( i i , j j ) = mean( data ( : , 2 ) ) ;v ( i i , j j ) = mean( data ( : , 3 ) ) ;w( i i , j j ) = mean( data ( : , 4 ) ) ;endendQ = sum( u ) ;C.29 profileQ.m%P r o f i l e s XY of the cross sec t ion wi th the ADVQdata = s t r u c t ( ' x ' ,{} , ' y ' ,{} , ' u ' ,{} ) ;sca le = .278 / 35000; %m/ stepmax width = 42000; %steps180max height = 60000; %stepsnum width samples = 4;num height samples = 9;f o l d e r = ' QPro f i l es / ' ;ExID = 7;ExTime = 0;%f i lename = [ fo l de r , ' E ' , num2str ( ExID ) ]%i f ( f i l e a t t r i b ( f i lename ) ==0)% mkdir ( f i lename ) ;%end%Send the ADV to home (0 ,0 )% pause ( travmove (5 ,0 ) ) ;% [ x , y ] = t ravgetpos ( )% pause ( travmove ( - x , 0 ) ) ;% [ x , y ] = t ravgetpos ( ) ;% pause ( travmove (0 , - y ) ) ;% [ x , y ] = t ravgetpos ( )Ve lMat r i x = [ ] ;x i = 0 ;y i = 0 ;t i cf p r i n t f ( t r avSe r i a lOb jec t , ' l ' ) ;se r i n = 'ReadyNOT ' ;wh i le ( ~ ( strcmp ( [ ' HomingComplete ' ,13 ,10 ] , se r i n ) ) )se r i n = fscan f ( t r a v S e r i a l O b j e c t ) ;end% i f ( toc > 2)% e r r o r ( ' Somethings up wi th t ravers ' )% endxstage = max width / ( num width samples - 1 ) ;s t r u c t I = 0 ;f o r i i = 1 : num height samplesf o r j j = 1 : num width samplespause ( 2 ) ;ADVgetdata ;Vdata = [ ] ;pause ( 0 . 5 ) ;f o r t =1:16ADVgetdata ;Vdata = [ Vdata ; ADVdataBuffer ] ;pause ( . 5 )ends t r u c t I = s t r u c t I + 1 ;Qdata ( s t r u c t I ) . x = x i ;181Qdata ( s t r u c t I ) . y = y i ;Qdata ( s t r u c t I ) . u = Vdata ;% f i lename = [ fo l de r , ' E ' , num2str ( ExID ) , ' / E ' , num2str ( ExID ) , ' - T ' , num2str (ExTime ) , ' -X ' , num2str ( x i ) , ' -Y ' , num2str ( y i ) , ' . t x t ' ]% save ( f i lename , ' data ' , ' - a s c i i ' , ' - tabs ' ) ;i f (mod( i i , 2 ) ==0)xstage = - round ( max width / ( num width samples - 1 ) ) ;e lsexstage = round ( max width / ( num width samples - 1 ) ) ;endx iy ii f ( j j < num width samples )x i = x i +xstage ;pause ( travmove ( xstage , 0 ) ) ;endendi f ( i i < num height samples )y i = y i + round ( max height / ( num height samples - 1 ) ) ;pause ( travmove (0 , round ( max height / num height samples ) ) ) ;endend%f i lename = [ fo l de r , ' E ' , num2str ( ExID ) , ' / E ' , num2str ( ExID ) , ' - T ' , num2str ( ExTime ) , ' .mat ' ]%save ( f i lename , ' scale ' , ' max width ' , ' max height ' , ' num width samples ' , 'num height samples ' ) ;C.30 SetupExperiment.mex = s t r u c t ( ' ExperimentNumber ' ,{} , ' Date ' ,{} , ' Note ' ,{} , 'Code ' ,{} ) ;sams = s t r u c t ( ' SampleNumber ' ,{} , ' ExperimentNumber ' ,{} ' Time ' ,{} , ' Note ' ,{} , 'Raw WL ',{} , 'Raw ADV ' ,{} , 'Q ' ,{} , ' dHt ' ,{} , ' ShutterPos ' ,{} , ' Weir1Pos ' ,{} , . . .' Weir2Pos ' ,{} , ' dHchanSet ' ,{} , ' Prc WL ' ,{} , ' Prc ADV ' ,{} ) ;%####### PARAMETERS #############shu t te rSe t = [0 ,3E6 , 3 . 5 E6,4E6 ] ;%####### SYSTEM INITIALIZATION ##########i f ( e x i s t ( ' we i r1Ser ia lOb jec t ' )&e x i s t ( ' s h u t t e r S e r i a l O b j e c t ' ) & e x i s t ( 't r a v S e r i a l O b j e c t ' ) )f p r i n t f ( ' System Setup Appears Okay . . . \ n ' ) ;e lsef p r i n t f ( ' I n i t i a l i z i n g System . . . ' ) ;s y s t e m i n i tsam i = 1 ;f p r i n t f ( ' I n i t i a l i z i n g System . . . [DONE] ' ) ;end182%########## RUN EXPERIMENT ###########f p r i n t f ( ' S t a r t i n g Experiment . . . \n ' ) ;f o r i i = 1 : leng th ( shu t te rSe t )endC.31 shutterinit.mf p r i n t f ( ' I n i t i a l i z i n g and homing s h u t t e r . . . ' )s h u t t e r S e r i a l O b j e c t = s e r i a l ( 'COM5 ' , ' BaudRate ' , 9600) ;fopen ( s h u t t e r S e r i a l O b j e c t ) ;f p r i n t f ( shu t t e rSe r i a lOb jec t , 13 )f p r i n t f ( shu t t e rSe r i a lOb jec t , ' mr 5000000 ' ) %Sends the s h u t t e r to f u l l open homep o s i t i o nf p r i n t f ( shu t t e rSe r i a lOb jec t , 13 ) ;pause (11) ;f p r i n t f ( shu t t e rSe r i a lOb jec t , 13 ) ;f p r i n t f ( shu t t e rSe r i a lOb jec t , ' p=0 ' ) %Sends the s h u t t e r to f u l l open home p o s i t i o nf p r i n t f ( shu t t e rSe r i a lOb jec t , 13 ) ;f p r i n t f ( ' [DONE]\n ' )% Vm = 400000 , iC.32 shuttermove.mf u n c t i o n shuttermove ( p o s i t i o n )g loba l s h u t t e r S e r i a l O b j e c t ;f p r i n t f ( shu t t e rSe r i a lOb jec t , 13 )f p r i n t f ( shu t t e rSe r i a lOb jec t , 'ma ' )f p r i n t f ( shu t t e rSe r i a lOb jec t , i n t 2 s t r ( - p o s i t i o n ) )f p r i n t f ( shu t t e rSe r i a lOb jec t , 13 )endC.33 smooth.mf u n c t i o n out = smooth ( in , width )i f (mod( width , 2 ) ==0)width = width - 1 ;endi f ( s i ze ( in , 2 )>1)183f p r i n t f ( 'WARNING: Smooth should only be used wi th a 1D ar ray .\n ' ) ;endband = ( width - 1 ) / 2 ;out = [ ] ;%Expanding band f o r the f i r s t measurementsf o r i = 1 : bandout = [ out ; mean ( 1 : ( 2 * i - 1 ) ) ] ;endf o r i =(band+1) : ( leng th ( i n ) -band )%[ i , ( i - band ) : ( i +band ) ]out = [ out ; mean( i n ( ( i - band ) : ( i +band ) ) ) ] ;end%Col laps ing band f o r the l a s t measurementsf o r i = ( leng th ( i n ) -band+1) : leng th ( i n )out = [ out ; mean( i n ( ( i - ( 2 * ( leng th ( i n ) - i ) ) ) : l eng th ( i n ) ) ) ] ;endC.34 startDAQ.mf u n c t i o n startDAQ ( )g loba l daqobject ;i f ( strcmp ( get ( daqobject , ' Running ' ) , 'On ' ) )wr ( 'DAQ Already running . ' ) ;s top ( daqobject ) ;endnewDAQdata ( ) ; %Clears b u f f e rs t a r t ( daqobject ) ;wr ( 'DAQ Star ted wi th c l ea r b u f f e r . ' ) ;C.35 systeminit.m%SYSTEM INITIALIZATIONi f ( e x i s t ( ' Sys temIn i t i a lT ime ' ) )wr ( 'NOTE: System prev ious l y i n i t i a l i z e d . ' ) ;e lseSys temIn i t i a lT ime = c lock ;d e f i n e g l o b a l ;UScal ib = load ( ' USCal ib ra t ion . t x t ' )%Weir 1w e i r 1 i n i t184%Weir 2w e i r 2 i n i t%Weir 2%Traverset r a v i n i t%Shut te rs h u t t e r i n i twr ( ' Shut te r p o s i t i o n must be set to zero when f u l l y open ' )%ADVADVini t%DAQd a q i n i t%Messageswr ( 'Move the wei r and update weir1zero to r e f l e c t new p o s i t i o n ' )endC.36 travgetpos.mf u n c t i o n [ x , y ] = t ravgetpos (sLDV)g loba l t r a v S e r i a l O b j e c t ;i f ( get ( t r avSe r i a lOb jec t , ' By tesAva i lab le ' ) ==0)f p r i n t f ( 'No p o s i t i o n data a v a i l b l e !\n ' ) ;x = NaN;y=NaN;r e t u r n ;endwhi le ( get ( t r avSe r i a lOb jec t , ' By tesAva i lab le ' )>0)se r i n = fscan f ( t r a v S e r i a l O b j e c t ) ;endindex = f i n d ( se r i n == ' , ' ) ;x = s t r2doub le ( se r i n ( 1 : ( index ( 1 ) -1) ) ) ;y = s t r2doub le ( se r i n ( ( index ( 1 ) +1) : ( index ( 2 ) -1) ) ) ;checksum = st r2doub le ( se r i n ( ( index ( 2 ) +1) : leng th ( se r i n ) ) ) ;i f ((1000000 - x - y ) ~= checksum )f p r i n t f ( ' Communications Er ro r w i th P o s i t i o n e r !\n\n ' ) ;endend185C.37 travinit.m%setups up communication wi th LDV P o s i t i o n e r%f c l o s e (sLDV)f p r i n t f ( ' Homing Traverse . . . ' )t r a v S e r i a l O b j e c t = s e r i a l ( 'COM3 ' , ' BaudRate ' , 57600) ;fopen ( t r a v S e r i a l O b j e c t ) ;se r i n = 'ReadyNOT ' ;wh i le ( ~ ( strcmp ( [ ' Ready . . . ' ,13 ,10 ] , se r i n ) ) )se r i n = fscan f ( t r a v S e r i a l O b j e c t ) ;endf p r i n t f ( t r avSe r i a lOb jec t , ' l ' ) ;wh i le ( ~ ( strcmp ( [ ' HomingComplete ' ,13 ,10 ] , se r i n ) ) )se r i n = fscan f ( t r a v S e r i a l O b j e c t ) ;endf p r i n t f ( ' [DONE]\n ' )C.38 travmove.mf u n c t i o n t =travmove ( x , y )g loba l t r a v S e r i a l O b j e c t ;speed = 30000/14; %steps / secondi f ( ( abs ( x )>0)&&(abs ( y ) ) )' cannot move i n both the x and y d i r e c t i o n s 'endi f ( abs ( x )>0)command = [ ' x ' , i n t 2 s t r i n g ( x ) ] ;e lsecommand = [ ' y ' , i n t 2 s t r i n g ( y ) ] ;endf p r i n t f ( t r avSe r i a lOb jec t , command) ;t =max( abs ( [ x , y ] ) ) * 1 . 1 / speed ;endC.39 ultraClean.mf u n c t i o n out = u l t raC lean ( i n )%Takes the daq ar ray i n and f i x e s a l l the channels%c lea r NaN from setout = [ ] ;nanitems = f i n d ( isnan ( i n ) ) ;186i f ( l eng th ( nanitems>0) )i f ( nanitems ( 1 ) > 1)nanitems = [ 0 ; nanitems ] ;endf o r i =2: leng th ( nanitems )i f ( nanitems ( i )<l eng th ( i n ) )[ ( nanitems ( i - 1 ) +1) : ( nanitems ( i ) -1) ] ;out = [ out ; i n ( ( nanitems ( i - 1 ) +1) : ( nanitems ( i ) -1) , : ) ] ;e lseout = [ out ; i n ( ( nanitems ( i - 1 ) +1) : leng th ( i n ) , : ) ] ;break ;endendelseout= i n ;end% f o r i = 1 : leng th ( i n )% i f (max( isnan ( i n ( i , : ) ) ) ==0)% out =[ out ; i n ( i , : ) ] ;% end% end%%This needs to be f i x e d .d t = [ ] ;f o r i =2:50d t = out ( i , 1 ) - out ( i - 1 ,1 ) ;endsample f req = 1/ min ( d t ) ;%sample f req = leng th ( out ) /max( out ( : , 1 ) ) ;f o r sens = 2: s ize ( out , 2 )%f i g u r e ( )% p l o t ( i n ( : , sens ) ) ;%Clear the spikes o f f f o r l oos ing contac tbase l ine = f f t smoo th ( out ( : , sens ) , sample freq , 0 , . 2 ) ;basel ine2 = f f t smoo th ( out ( : , sens ) , sample freq , 0 , . 0 5 ) ;out2=out ;th resho ld = 0 . 2 ;sp ikes = f i n d ( out ( : , sens )>(base l ine+ th resho ld ) ) ;f o r t = 1 : leng th ( sp ikes )i f ( sp ikes ( t ) / l eng th ( out ) < 0.95)out ( sp ikes ( t ) , sens ) = basel ine2 ( spikes ( t ) ) ;endend%% subp lo t (3 ,1 ,1 ) ;% p l o t ( out2 ( : , sens ) - basel ine , ' r ' ) ;% hold on ;% p l o t ( out ( : , sens ) - base l ine ) ;187% subp lo t (3 ,1 ,2 ) ;% p l o t ( out2 ( : , sens ) , ' r ' ) ;% hold on% p l o t ( out ( : , sens ) , ' b ' ) ;% p l o t ( basel ine , ' g ' ) ;% p l o t ( basel ine2 , ' c ' ) ;%% % Low Pass f i l t e r% subp lo t (3 ,1 ,3 ) ;% out ( : , sens ) = f f t smoo th ( out ( : , sens ) , sample freq , 0 , 1 ) ;% p l o t ( out ( : , 1 ) , out ( : , sens ) , ' b ' ) ;endC.40 USdepth.mf u n c t i o n d = USdepth (V,SensNum)g loba l UScal ibi f ( ( ( SensNum <= 2) &(V > 4 .9 ) ) | ( ( SensNum == 1) &(V > 6 .5 ) ) )d = nan ;d isp ( 'NAN f o r depth ' ) ;e lse%% USE CALIBRATION TO GET DISTANCE TO BOTTOM OF FLUMEset= f i n d ( UScal ib ( : , 1 ) ==SensNum) ;i f ( l eng th ( set ) ~=2)e r r o r ( ' c a l i b r a t i o n s ize wrong ' ) ;endc = UScal ib ( set , 2 : 3 ) ;m = ( c (2 ,2 ) - c (1 ,2 ) ) / ( c (2 ,1 ) - c (1 ,1 ) ) ;b = c (1 ,2 ) - m* c (1 ,1 ) ;d = V . *m+b ;%APPLY SURVEY POINTSd = d+( UScal ib ( ( SensNum*2) ,4 ) - min ( UScal ib ( : , 4 ) ) ) /100 ;endC.41 usZero.m%usZero%Run f o r 2 minutes wi th no f low . Ca lcu la tes us1zero , us2zero , us3zero188%Clear AI Bu f fe ri f ( get ( daqobject , ' SamplesAvai lable ' )>0)getdata ( daqobject , get ( daqobject , ' SamplesAvai lable ' ) ) ;ends t a r t ( daqobject )t i c ;zerodata = [ ] ;wh i le ( toc<120)i f ( get ( daqobject , ' SamplesAvai lable ' )>0)zerodata = [ zerodata ; getdata ( daqobject , get ( daqobject , ' SamplesAvai lable ' ) )] ;endi f (mod( toc , 5 )<0.1)f p r i n t f ( i n t 2 s t r ( toc ) )f p r i n t f ( '\n ' )endpause ( 0 . 0 1 )endstop ( daqobject )i f ( get ( daqobject , ' SamplesAvai lable ' )>0)zerodata = [ zerodata ; getdata ( daqobject , get ( daqobject , ' SamplesAvai lable ' ) ) ] ;endp l o t ( zerodata ( : , 1 ) ) ;hold on ;p l o t ( zerodata ( : , 2 ) , ' r ' ) ;p l o t ( zerodata ( : , 3 ) , ' g ' )us1zero = mean( zerodata ( : , 1 ) ) ;us2zero = mean( zerodata ( : , 2 ) ) ;us3zero = mean( zerodata ( : , 3 ) ) ;C.42 weir1init.mwe i r1Ser ia lF lag = 0;weir1zero =0;w e i r 1 s t a r t s e r i a lC.43 weir1move.mf u n c t i o n Move2Pos ( pos )g loba l we i r1Ser ia lOb jec t ;g loba l weir1zero ;% ' moving '% f p r i n t f ( num2str ( round ( pos ) ) ) ;f p r i n t f ( we i r1Ser ia lOb jec t , num2str ( round ( pos ) ) ) ;189endC.44 weir1startserial.m%opens communication wi th Arduino a t Weir 1 .i f ( we i r1Ser ia lF lag ~= 1)we i r1Ser ia lF lag = 1;we i r1Ser ia lOb jec t = s e r i a l ( 'COM4 ' , ' BaudRate ' , 9600) ;fopen ( we i r1Ser ia lOb jec t ) ;e lse' Cannot s t a r t - wei r 1 a l read i n i t i a l i z e d . ( Flag Set ) 'end%%% i f ( we i r1Ser ia lF lag ~= 1)% we i r1Ser ia lF lag = 1;% we i r1Ser ia lOb jec t = s e r i a l ( 'COM4' , ' BaudRate ' , 57600) ;% fopen ( we i r1Ser ia lOb jec t ) ;% else% ' Cannot s t a r t - wei r 1 a l read i n i t i a l i z e d . ( Flag Set ) '% endC.45 weir2init.mwe i r2Ser ia lF lag = 0;weir2zero =0;w e i r 2 s t a r t s e r i a lC.46 weir2move.mf u n c t i o n Move2Pos ( pos )g loba l we i r2Ser ia lOb jec t ;g loba l weir2zero ;f p r i n t f ( we i r2Ser ia lOb jec t , num2str ( round ( pos ) ) ) ;end190C.47 weir2startserial.m%opens communication wi th Arduino a t Weir 1 .i f ( we i r2Ser ia lF lag ~= 1)we i r2Ser ia lF lag = 1;we i r2Ser ia lOb jec t = s e r i a l ( 'COM6 ' , ' BaudRate ' , 9600) ;fopen ( we i r2Ser ia lOb jec t ) ;e lse' Cannot s t a r t - wei r 1 a l read i n i t i a l i z e d . ( Flag Set ) 'end191Appendix DMatlab Code for HydrokineticModelD.1 CalcPower.mf u n c t i o n ResultsMat = CalcPower ( l oca t i ons , Areas )g loba l channe l mat r i x ;g loba l ResultsMat ;g loba l pTurbA ;f o r i = 1 : leng th ( l o c a t i o n s )channe l mat r i x ( l o c a t i o n s ( i ) , pTurbA ) = Areas ( i ) ;endsub sub ( ) ;ResultsMat ( : , : , s i ze ( ResultsMat , 3 ) +1) = channe l mat r i x ;endD.2 QReduction.mchanne l mat r i x ( : , pTurbA ) =0;%[ PowerSum , channe l mat r i x ] = sub sub ( ) ;t i c192Aset = 200;QSet = [ 9 0 . 9 , 1 0 0 ] ;ResultsMat = [ ] ;f o r j = 1 : leng th ( QSet )channe l mat r i x ( : , pQ) = QSet ( j ) ;f o r i = 1 : leng th ( Aset )i f ( Aset ( i ) <= channe l mat r i x (3300 , pArea ) )CalcPower (3300 , Aset ( i ) ) ;e l s e i f ( Aset ( i ) <= ( channe l mat r i x (3300 , pArea ) +channe l mat r i x (3100 , pArea ) ))CalcPower ( [3300 ,3100 ] , [ channe l mat r i x (3300 , pArea ) , Aset ( i ) -channe l mat r i x (3300 , pArea ) ] ) ;e l s e i f ( Aset ( i ) <= ( channe l mat r i x (3300 , pArea ) +channe l mat r i x (3100 , pArea ) +channe l mat r i x (2900 , pArea ) ) ) ;CalcPower ( [3300 ,3100 ,2900] , [ channe l mat r i x (3300 , pArea ) , channe l mat r i x(3100 , pArea ) , Aset ( i ) - channe l mat r i x (3300 , pArea ) - channe l mat r i x(3100 , pArea ) ] ) ;e lseCalcPower ( [3300 ,3100 ,2900 ,2700] , [ channe l mat r i x (3300 , pArea ) ,channe l mat r i x (3100 , pArea ) , channe l mat r i x (2900 , pArea ) , Aset ( i ) -channe l mat r i x (3300 , pArea ) - channe l mat r i x (3100 , pArea ) -channe l mat r i x (2900 , pArea ) ] ) ;endendendResultsMat = ResultsMat ( : , : , 2 : s i ze ( ResultsMat , 3 ) ) ;tocD.3 get properties.mf u n c t i o n [A, Rh, Wp] = s e t p r o p e r t i e s ( index , channe l mat r i x )[ channe l mat r i x ( index , pArea ) , channe l mat r i x ( index , pRh) , channe l mat r i x ( index,pWp) ] = s e c t i o n p r o p e r t i e s ( channe l mat r i x ( index , pc lass ) , channe l mat r i x (index , pdepth ) ) ;end193D.4 smoothClose.mf u n c t i o n i n =smoothClose ( in , d i s t )f o r i = 1 : leng th ( i n )myset = f i n d ( abs ( i n ( : , 1 ) - i n ( i , 1 ) )<d i s t ) ;i n ( i , 2 ) = mean( i n ( myset , 2 ) ) ;endD.5 CalculateLength.mchanLength = zeros ( leng th ( r ) ,1 ) ;f o r i = 2 : leng th ( r )chanLength ( i ) = chanLength ( i - 1 ) + s q r t ( ( r ( i , 1 ) - r ( i - 1 ,1 ) ) ˆ2 + ( r ( i , 2 ) - r ( i - 1 ,2 )) ˆ 2 ) ;endo f f s e t =( r ( : , 4 ) - chanLength ) ;out = [ ] ;f o r i = 1 : leng th ( o f f s e t )i f ( ~ isnan ( o f f s e t ( i ) ) )out = o f f s e t ( i ) ;endendmean( out )D.6 RemoveData.mf u n c t i o n Out=RemoveData ( In , Thresh )load k i l l d a t a . t x tf o r i = 1 : leng th ( In )myset = f i n d ( abs ( In ( i , 1 ) - k i l l d a t a ( : , 1 ) ) < Thresh ) ;i f ( ~ isempty ( myset ) )myset2 = f i n d ( abs ( In ( i , 2 ) - k i l l d a t a ( myset , 2 ) ) < Thresh ) ;i f ( ~ isempty ( myset2 ) )In ( i , : ) = NaN;endendendOut = [ ] ;194f o r i = 1 : leng th ( In ) ;i f ( ~ isnan ( In ( i , 1 ) ) )Out = [ Out ; In ( i , : ) ] ;endendD.7 range of TurbA.mchanne l mat r i x ( : , pTurbA ) =0;[ PowerSum , channe l mat r i x ] = sub sub ( ) ;t i c%Aset = 0 .0 :0 .0005 :0 .005 ;%Aset = [ 0 : 5 : 5 0 0 ] ;Aset = [ 5 0 2 : 2 : 8 0 0 ] ;%ResultsMat = [ ] ;f o r i = 1 : leng th ( Aset )i / l eng th ( Aset )i f ( Aset ( i ) <= channe l mat r i x (3300 , pArea ) )CalcPower (3300 , Aset ( i ) ) ;e l s e i f ( Aset ( i ) <= ( channe l mat r i x (3300 , pArea ) +channe l mat r i x (3100 , pArea ) ) )CalcPower ( [3300 ,3100 ] , [ channe l mat r i x (3300 , pArea ) , Aset ( i ) - channe l mat r i x(3300 , pArea ) ] ) ;e l s e i f ( Aset ( i ) <= ( channe l mat r i x (3300 , pArea ) +channe l mat r i x (3100 , pArea ) +channe l mat r i x (2900 , pArea ) ) ) ;CalcPower ( [3300 ,3100 ,2900] , [ channe l mat r i x (3300 , pArea ) , channe l mat r i x(3100 , pArea ) , Aset ( i ) - channe l mat r i x (3300 , pArea ) - channe l mat r i x (3100 ,pArea ) ] ) ;e lseCalcPower ( [3300 ,3100 ,2900 ,2700] , [ channe l mat r i x (3300 , pArea ) ,channe l mat r i x (3100 , pArea ) , channe l mat r i x (2900 , pArea ) , Aset ( i ) -channe l mat r i x (3300 , pArea ) - channe l mat r i x (3100 , pArea ) - channe l mat r i x(2900 , pArea ) ] ) ;endend%ResultsMat = ResultsMat ( : , : , 2 : s i ze ( ResultsMat , 3 ) ) ;ProcessResultsMattocbeep ;pause ( 0 . 5 )beep ;pause ( 0 . 5 )beep ;pause ( 1 )195beep ;D.8 solve flow.mf u n c t i o n channe l mat r i x = so l ve f l ow ( channe l mat r i x )number ce l ls = max( s ize ( channe l mat r i x ) ) ;%St ruc tu re o f channe l mat r i x%Class , kf , k t , Depth , Ve loc i t y , Power , Lengtht o t a l h e a d d i f f e r e n c e = channe l mat r i x ( number cel ls , pdepth ) - channe l mat r i x(1 , pdepth )[ A1 , Rh1 , Wp1] = g e t p r o p e r t i e s (1 , channe l mat r i x )[ A2 , Rh2 , Wp2] = g e t p r o p e r t i e s ( number cel ls , channe l mat r i x )AveA = (A1+A2) / 2 ;AveRh = (Rh1+Rh2) / 2 ;channe l mat r i x = reca lcu la te channe l params ( channe l mat r i x ) ;%Clacu la te i n i t i a l guess a t f low ra te based on mean p r o p e r t i e sQ = s q r t ( t o t a l h e a d d i f f e r e n c e * 2 * g * mean( channe l mat r i x ( : , pRh) ) * (mean( channe l mat r i x ( : , pArea ) ) ) ˆ2 / sum( channe l mat r i x ( : , p length ) ) ) ;f i r s t l o o p = 1;f o r counter = 1 : number cel ls -1i =number cel ls - counter +1;channe l mat r i x ( i , pve l ) = Q/ channe l mat r i x ( i , pArea ) ;delh = channe l mat r i x ( i , pk f ) * ( channe l mat r i x ( i , p length ) /channe l mat r i x ( i , pRh) ) * channe l mat r i x ( i , pve l ) ˆ2 / (2 * g ) ;channe l mat r i x ( i -1 , pTotalHead ) = channe l mat r i x ( i -1 , pTotalHead ) +delh ;channe l mat r i x ( i - 1 , : ) = r e c a l c c e l l ( channe l mat r i x ( i - 1 , : ) ) ; %r e c a l c u l a t e area , rh , e tc .endend196D.9 InterpolateDepths.mload surveyDepths . t x t ;surveyDepths ( : , 1 ) = surveyDepths ( : , 1 ) - surveyDepths (1 ,1 ) ;number of po in ts = 3562;e l ev a t i on = zeros (1 , number of po in ts ) ;e l ev a t i on ( 1 ) = surveyDepths ( 1 ) ;marker = 2 ;f o r i = 2 : number of po in tsi f ( surveyDepths ( marker , 1 ) == i )%I f t h i s c e l l i s on a survey po in t , use t h a t .e l ev a t i on ( i ) = surveyDepths ( marker , 2 ) ;marker = marker +1;e lse%L i n e a r l y I n t e r p o l a t e Between Poin tsslope = ( surveyDepths ( marker , 2 ) - surveyDepths ( marker -1 ,2 ) ) / ( surveyDepths( marker , 1 ) - surveyDepths ( marker -1 ,1 ) ) ;e l ev a t i on ( i ) = ( i - surveyDepths ( marker -1 ,1 ) ) * s lope + surveyDepths ( marker-1 ,2 ) ;endendp l o t ( surveyDepths ( 2 : leng th ( surveyDepths ) ,1 ) , surveyDepths ( 2 : leng th ( surveyDepths ),2 ) , ' r ' ) ;hold on ;%p l o t ( 1 : number of po ints , e leva t ion , ' b ' ) ;D.10 SensativityPlot.mf i g u r e ( )se t ( gcf , ' Un i ts ' , ' Inches ' )se t ( gcf , ' Pos i t i on ' , [ 0 , 0 , 6 , 4 . 5 ] )%subp lo t (2 ,1 ,1 ) ;se t ( gca , ' FontName ' , ' Pa la t i no ' )hold a l l%Flowrange = [ 2 , 1 , 3 , 4 ] ;index = 2;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' k-> ' )%Turbine C o e f f i c i e n trange = [ 7 , 1 , 8 ] ;index = 3;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' g - - . ' )197%Depthrange = [ 5 , 1 , 6 ] ;index = 4;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' b : o ' )legend ( ' Flow ' , ' Turbine Drag Rat io ' , ' Measured Depth ' )%subp lo t (2 ,1 ,2 ) ;se t ( gca , ' FontName ' , ' Pa la t i no ' )hold a l l%Flowrange = [ 2 , 1 , 3 , 4 ] ;index = 2;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 8 ) , ' k-> ' )%Turbine C o e f f i c i e n trange = [ 7 , 1 , 8 ] ;index = 3;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 8 ) , ' g - - . ' )%Depthrange = [ 5 , 1 , 6 ] ;index = 4;p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 8 ) , ' b : o ' )%% %Flow% range = [ 2 , 1 , 3 , 4 ] ;% index = 2;% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' k . ' )% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' ko ' )%%%% %Turbine C o e f f i c i e n t% range = [ 7 , 1 , 8 ] ;% index = 3;% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' k . ' )% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' ko ' )%% %Depth% range = [ 5 , 1 , 6 ] ;% index = 4;% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' k . ' )% p l o t ( Sens ( range , index ) . / Sens (1 , index ) , Sens ( range , 7 ) , ' ko ' )198D.11 range of kt.mf p r i n t f ( 'No Turbines ' ) ;channe l mat r i x ( : , pk t ) =0;[ PowerSum , channe l mat r i x ] = sub sub ( ) ;[ PowerSum , channe l mat r i x ] = sub sub ( ) ;Performance = [0 , channe l mat r i x (2 , pdepth ) ] ;f i g u r e ( ) ;p l o t ( channe l mat r i x ( : , pe leva t i on ) , ' k ' ) ;hold on ;p l o t ( channe l mat r i x ( : , pdepth ) +channe l mat r i x ( : , pe leva t i on ) , ' b ' )k tSet = 0 .0005 :0 .0005 :0 .01 ;f o r i = 1 : leng th ( k tSet )% se tup mat r i xt h i s k t = k tSet ( i ) /10channe l mat r i x (3100:3200 , pk t ) = t h i s k t ;% p l o t ( channe l mat r i x ( : , pk t ) , ' r ' )[ TurbineHead , channe l mat r i x ] = sub sub ( ) ;Performance = [ Performance ; TurbineHead , channe l mat r i x (2 , pdepth ) ] ;p l o t ( channe l mat r i x ( : , pdepth ) +channe l mat r i x ( : , pe leva t i on ) , ' r ' )%p l o t ( channe l mat r i x ( : , pk t ) , ' . ' )endf i g u r e ( )p l o t ( Performance ( : , 1 ) *998*9.81*110* .3 , ' g ' )hold on ;p l o t ( ( - Performance ( : , 2 ) +Performance (1 ,2 ) ) *998*9.81*110* .6 , ' r ' )p l o t ( ( Performance ( : , 1 ) * . 3 - Performance ( : , 2 ) * .6+ Performance (1 ,2 ) * . 6 ) *998*9.81*110 , 'b ' )legend ( ' Hyd rok ine t i c Generat ion ' , ' Loss o f Large Hydro Generat ion ' , ' Net Generat ion' ) ;x l a b e l ( ' Number o f t u rb i nes (Low Resistance Range ) ' ) ;y l a b e l ( ' Power (w) ' ) ;D.12 sub sub.mf u n c t i o n [ PowerSum , channe l mat r i x ]= sub sub ( )%Takes the channel mat r i x and does backf low c a l c u l a t i o n s . E n t i r e sec t ion ds%c o n t r o l l e d .199g loba l channe l mat r i x ;g loba l rho ;g loba l g ;g loba l pc lass ;g loba l pk f ;g loba l pk t ;g loba l pdepth ;g loba l pve l ;g loba l ppower ;g loba l p length ;g loba l pRh ;g loba l pWp;g loba l pArea ;g loba l pTotalHead ;g loba l pe leva t i on ;g loba l pQ;g loba l ma t r i x w id th ;g loba l pFr ;g loba l pX ;g loba l pTurbA ;g loba l f i l e n o m i n a l d e p t h s ;g loba l w e t t e d p a r i m e t e r s f r o m f i l e ;g loba l a r e a s f r o m f i l e ;indexSet = s o r t ( 2 : leng th ( channe l mat r i x ) , ' descend ' ) ;PowerSum = 0;f o r i = indexSet%Using t o t a l head , determine v e l o c i t ychanne l mat r i x ( i , : ) = r e c a l c c e l l ( channe l mat r i x ( i , : ) ) ;i f ( channe l mat r i x ( i , pTurbA )>0)Power = 0 .5*0 .25* rho * channe l mat r i x ( i , pve l ) ˆ3 * channe l mat r i x ( i , pTurbA ) ;channe l mat r i x ( i , ppower ) = Power ;delh = Power / (0 .5 * rho * g * channe l mat r i x ( i ,pQ) ) ;PowerSum = PowerSum + Power ;e lsedelh = 0;enddelh = delh + channe l mat r i x ( i , pk f ) * ( channe l mat r i x ( i , p length ) / channe l mat r i x( i , pRh) ) * channe l mat r i x ( i , pve l ) ˆ2 / (2 * g ) ;channe l mat r i x ( i -1 , pTotalHead ) = channe l mat r i x ( i , pTotalHead ) + delh ;endPowerSumchanne l mat r i x ( 1 , : ) = r e c a l c c e l l ( channe l mat r i x ( 1 , : ) ) ;D.13 LoadDataFiles.m%Loads data from f i l e s . Filenames are set i n ModelSetup .m200f i l e n o m i n a l d e p t h s = load ( DataFi leDepths ) ;w e t t e d p a r i m e t e r s f r o m f i l e = load ( WPFile ) ;a r e a s f r o m f i l e = load ( AreaF i le ) ;bed e leva t ion = load ( BedEleva t ionF i le ) ;D.14 SolveKf.mf u n c t i o n k f = SolveKf ( depth )%Solves f o r the f r i c t i o n c o e f f i c i e n t to achieve a uni form depth Upstream%and Downstream .f o r x=1g loba l rho ;g loba l g ;g loba l pc lass ;g loba l pk f ;g loba l pk t ;g loba l pdepth ;g loba l pve l ;g loba l ppower ;g loba l p length ;g loba l pRh ;g loba l pWp;g loba l pArea ;g loba l pTotalHead ;g loba l pe leva t i on ;g loba l pQ;g loba l ma t r i x w id th ;g loba l pFr ;g loba l pX ;g loba l f i l e n o m i n a l d e p t h s ;g loba l w e t t e d p a r i m e t e r s f r o m f i l e ;g loba l a r e a s f r o m f i l e ;g loba l channe l mat r i x ;end %Global Var iab les . In f o r loop to hide .channe l mat r i x ( : , pdepth ) = depth ;channe l mat r i x ( : , pk t ) = 0 ;k f = 0 .02 ;e r r o r = 1E10 ; %Something la rge to ensure a t l e a s t one loop .wh i le ( abs ( e r r o r ) > 0.05)channe l mat r i x ( : , pdepth ) = depth ;channe l mat r i x ( : , pk f ) = k f ;sub sub ( ) ;e r r o r = channe l mat r i x (1 , pdepth ) - depthk f = k f * ( 1 - ( 2 * e r r o r / depth ) )201endD.15 recalc cell.mf u n c t i o n c e l l d a t a = r e c a l c c e l l ( c e l l d a t a )pc lass =1;pk f =2;pk t =3;pdepth =4;pve l =5;ppower =6;p length = 7;pRh = 8;pWp =9;pArea =10;pTotalHead =11;pe leva t i on = 12;pQ=13;pFr = 14;pX = 15;ma t r i x w id th = 13;g=9.81;%Changes the depth , v e l o c i t y to match the t o t a l head .e r r o r = 1E10 ;% A la rge number to ensure a loop .wh i le ( e r r o r > 0.0001)c e l l d a t a ;c e l l d a t a ( pdepth ) = c e l l d a t a ( pTotalHead ) - c e l l d a t a ( pve l ) ˆ 2 / ( 2 * g ) -c e l l d a t a ( pe leva t i on ) ;i f ( c e l l d a t a ( pdepth ) < 0)c e l l d a t a ( pdepth ) = 0 . 5 ;f p r i n t f ( ' E r ro r - negat ive depth\n ' ) ;end[A, Rh, Wp] = s e c t i o n p r o p e r t i e s ( c e l l d a t a ( pc lass ) , c e l l d a t a ( pdepth ) ,c e l l d a t a (pX) ) ;c e l l d a t a ( pArea ) = A;c e l l d a t a (pRh) = Rh ;c e l l d a t a (pWp) = Wp;c e l l d a t a ( pve l ) = c e l l d a t a (pQ) /A ;e r r o r = ( c e l l d a t a ( pdepth ) + c e l l d a t a ( pve l ) ˆ 2 / ( 2 * g ) + c e l l d a t a (pe leva t i on ) ) - c e l l d a t a ( pTotalHead ) ;c e l l d a t a ( pdepth ) = c e l l d a t a ( pdepth ) + e r r o r ;c e l l d a t a ( pFr ) = c e l l d a t a ( pve l ) / s q r t ( c e l l d a t a ( pdepth ) *g ) ;endend202D.16 ModelSetup.mDataFi leDepths = ' SeatonDataf i leDepths . t x t ' ;WPFile = 'SeatonWP . t x t ' ;AreaF i le = ' SeatonArea . t x t ' ;BedEleva t ionF i le = ' SeatonBedElevation . t x t ' ;%NominalDepth = 6; %metersD.17 SolveKf real.mf u n c t i o n Compares= So lveK f rea l ( )%Solves f o r the f r i c t i o n c o e f f i c i e n t to achieve a uni form depth Upstream%and Downstream .f o r x=1g loba l rho ;g loba l g ;g loba l pc lass ;g loba l pk f ;g loba l pk t ;g loba l pdepth ;g loba l pve l ;g loba l ppower ;g loba l p length ;g loba l pRh ;g loba l pWp;g loba l pArea ;g loba l pTotalHead ;g loba l pe leva t i on ;g loba l pQ;g loba l ma t r i x w id th ;g loba l pFr ;g loba l pX ;g loba l f i l e n o m i n a l d e p t h s ;g loba l w e t t e d p a r i m e t e r s f r o m f i l e ;g loba l a r e a s f r o m f i l e ;g loba l channe l mat r i x ;end %Global Var iab les . In f o r loop to hide .g loba l NominalDepth ;Compares = [ ] ;f o r NominalDepth = [ 5 . 4 , 5 . 4 2 5 ]%[ 5 . 95 , 5 . 96 , 5 . 97 ]se tup mat r i x ;depth=NominalDepth ;channe l mat r i x ( : , pk t ) = 0 ;203k f = 0.00286;e r r o r = 1E10 ; %Something la rge to ensure a t l e a s t one loop .%whi le ( abs (nanmean ( e r r o r ) ) > 0.05)%f i g u r e ( ) ;%p l o t ( channe l mat r i x ( : , pMeasDepth ) , ' r ' ) ;%hold on ;f o r k f = 0.0020:0.0001:0.0025channe l mat r i x ( : , pdepth ) = depth ;channe l mat r i x ( : , pk f ) = k f ;sub sub ( ) ;e r r o r = ( channe l mat r i x ( : , pdepth ) - channe l mat r i x ( : , pMeasDepth ) * 0 . 9 ) ;%p l o t ( channe l mat r i x ( : , pdepth ) +channe l mat r i x ( : , pe leva t i on ) ) ;Compares = [ Compares ; NominalDepth , kf , nansum ( e r r o r . ˆ 2 ) ] ;%k f = k f * ( 1 - ( 2 * nanmean ( e r r o r ) / depth ) )%k f%cor rcoe f ( channe l mat r i x (515:3470 , pdepth ) , channe l mat r i x (515:3470 ,p%MeasDepth ) )endendD.18 recalculate channel params.mf u n c t i o n channe l mat r i x = reca lcu la te channe l params ( channe l mat r i x )f o r index =1:max( s ize ( channe l mat r i x ) ) ;channe l mat r i x ( index , : ) = r e c a l c c e l l ( channe l mat r i x ( index , : ) ) ;endendD.19 untitled6.mPower = [ ] ;HeadChange = [ ] ;TurbArea = [ ] ;f o r i = 1 : s ize ( ResultsMat , 3 )Power = [ Power ; sum( ResultsMat ( : , pTurbA , i ) . * ResultsMat ( : , pvel , i ) . ˆ 3 * 0 . 5 * . 2 5 *rho ) ] ;HeadChange = [ HeadChange ; ResultsMat (1 , pTurbA , i ) - ResultsMat (1 , pTurbA , 1 ) ]TurbArea = [ TurbArea = sum( ResultsMat ( : , pTurbA , i ) ] ;204D.20 MonteCarloTest.m% Example Monte Car lo S imu la t ion i n Matlab% Funct ion : y = x2 ˆ 2 / x1%% Generate n samples from a normal d i s t r i b u t i o n% r = ( randn ( n , 1 ) * sd ) + mu% mu : mean% sd : standard d e v i a t i o n%% Generate n samples from a uni form d i s t r i b u t i o n% r = a + rand ( n , 1 ) * ( b - a )% a : minimum% b : maximumn = 100000; % The number o f f u n c t i o n eva lua t ions% - - - Generate vec to rs o f random inpu ts% x1 ~ Normal d i s t r i b u t i o n N(mean=100 ,sd=5)% x2 ~ Uniform d i s t r i b u t i o n U( a=5 ,b=15)x1 = ( randn ( n , 1 ) * 5 ) + 100;x2 = 5 + rand ( n , 1 ) * ( 15 - 5 ) ;% - - - Run the s imu la t i on% Note the use of element - wise m u l t i p l i c a t i o ny = x2 . ˆ 2 . / x1 ;% - - - Create a histogram of the r e s u l t s (50 bins )h i s t ( y ,50 ) ;% - - - Ca lcu la te summary s t a t i s t i c sy mean = mean( y )y s td = std ( y )y median = median ( y )D.21 define global.mg loba l rho ;g loba l g ;g loba l pc lass ;g loba l pk f ;g loba l pk t ;g loba l pdepth ;g loba l pve l ;g loba l ppower ;g loba l p length ;205g loba l pRh ;g loba l pWp;g loba l pArea ;g loba l pTotalHead ;g loba l pe leva t i on ;g loba l pQ;g loba l ma t r i x w id th ;g loba l pFr ;g loba l pX ;g loba l pMeasDepth ;g loba l pTurbA ;g loba l f i l e n o m i n a l d e p t h s ;g loba l w e t t e d p a r i m e t e r s f r o m f i l e ;g loba l a r e a s f r o m f i l e ;g loba l channe l mat r i x ;g loba l ResultsMatg=9.81;rho =1000;%Def ine i n d i c i e s f o r channe l mat r i xpc lass =1;pk f =2;pk t =3;pdepth =4;pve l =5;ppower =6;p length = 7;pRh = 8;pWp =9;pArea =10;pTotalHead =11;pe leva t i on = 12;pQ=13;pFr = 14;pX = 15;pMeasDepth = 16;pTurbA = 17;ma t r i x w id th = 17;D.22 section properties.mf u n c t i o n [A, Rh, Wp] = s e c t i o n p r o p e r t i e s ( c lass , depth , c e l l p o s )g loba l f i l e n o m i n a l d e p t h s ;g loba l w e t t e d p a r i m e t e r s f r o m f i l e ;g loba l a r e a s f r o m f i l e ;i f ( c lass == 1)%This i s a t r i a n g u l a r channel , 2m deep and 5m accross a t top .A = 0.5 * depth * depth * 5 / 2 ;Wp = s q r t ( depth ˆ2 + ( depth * 5 / 2) ˆ 2 ) ;206Rh = A/Wp;e l s e i f ( c lass == 2)%This uses the data loaded from f i l e s and l i n e a r l y i n t e r p o l a t e s%between po in t s .Wp = i n t e r p 1 ( f i l e nom ina l dep ths , w e t t e d p a r i m e t e r s f r o m f i l e ( ce l l pos , : ) ,depth ) ;A = i n t e r p 1 ( f i l e nom ina l dep ths , a r e a s f r o m f i l e ( ce l l pos , : ) , depth ) ;Rh = A/Wp;endendD.23 PlotSensitivity.mf i g u r e ( )Performance = Perf2 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' k ' )hold on ;Performance = Perf3 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' k - - ' )Performance = Perf8 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' k - . ' )Performance = Perf4 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' k : ' )Performance = Perf5 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' k< ' )Performance = Perf6 ;p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) ) *9.81*110 ,' kx ' )legend ( ' 4m ' , ' 4.4m ' , ' 5.0m ' , ' 5.4m ' , ' 6.0m ' , ' 6.4m ' )x l a b e l ( ' Number o f Turbines ( High Resistance Range ) ' ) ;y l a b e l ( ' Net Power (kW) ' ) ;% Performance = Perf2 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' k ' )% hold on ;% Performance = Perf3 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' k - - ' )% Performance = Perf8 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )207*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' k - . ' )% Performance = Perf4 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' k : ' )% Performance = Perf5 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' k< ')% Performance = Perf6 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 - ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 ) )*9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 ) * .85 - Performance (1 ,2 ) * . 8 5 )) , ' kx ' )% f i g u r e ( )%% Performance = Perf2 ;% p l o t ( ( Performance ( : , 1 ) * . 3 5 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' k ' )% hold on ;% Performance = Perf3 ;% p l o t ( ( Performance ( : , 1 ) * .355 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' k - - ' )% Performance = Perf8 ;% p l o t ( ( Performance ( : , 1 ) * .355 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' k - . ' )% Performance = Perf4 ;% p l o t ( ( Performance ( : , 1 ) * .355 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' k : ' )% Performance = Perf5 ;% p l o t ( ( Performance ( : , 1 ) * .355 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' k< ')% Performance = Perf6 ;% p l o t ( ( Performance ( : , 1 ) * .355 ) *9 .81*110 , ( Performance ( : , 1 ) * . 3 5 5 . / ( Performance ( : , 2 )* .85 - Performance (1 ,2 ) * . 8 5 ) ) , ' kx ' )%% legend ( ' 4 . 4m' , ' 5 . 0m' , ' 5 . 4m' , ' 6 . 0m' , ' 6 . 4m' )% x l a b e l ( ' Hyd rok ine t i c Generat ion (kW) ' ) ;% y l a b e l ( ' Rat io o f HK Generat ion to Dam Generat ion Loss ' ) ;D.24 deg2utm.mf u n c t i o n [ x , y , utmzone ] = deg2utm ( Lat , Lon )% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -% [ x , y , utmzone ] = deg2utm ( Lat , Lon )%% Desc r i p t i on : Funct ion to conver t l a t / lon vec to rs i n t o UTM coord ina tes (WGS84) .% Some code has been ex t rac ted from UTM.m f u n c t i o n by Gabr ie l Ruiz Mart inez .%% Inpu ts :% Lat : La t i t ude vec to r . Degrees . +ddd . ddddd WGS84% Lon : Longi tude vec to r . Degrees . +ddd . ddddd WGS84208%% Outputs :% x , y , utmzone . See example%% Example 1 :% Lat =[40.3154333; 46.283900; 37.577833; 28.645650; 38.855550; 25.061783] ;% Lon=[ -3.4857166; 7.8012333; -119.95525; -17.759533; -94.7990166;121.640266];% [ x , y , utmzone ] = deg2utm ( Lat , Lon ) ;% f p r i n t f ( '%7.0 f ' , x )% 458731 407653 239027 230253 343898 362850% f p r i n t f ( '%7.0 f ' , y )% 4462881 5126290 4163083 3171843 4302285 2772478% utmzone =% 30 T% 32 T% 11 S% 28 R% 15 S% 51 R%% Example 2 : I f you have Lat / Lon coord ina tes i n Degrees , Minutes and Seconds% LatDMS=[40 18 55.56; 46 17 2 . 0 4 ] ;% LonDMS=[ -3 29 8 .58 ; 7 48 4 . 4 4 ] ;% Lat=dms2deg ( mat2dms (LatDMS) ) ; %conver t i n t o degrees% Lon=dms2deg ( mat2dms (LonDMS) ) ; %conver t i n t o degrees% [ x , y , utmzone ] = deg2utm ( Lat , Lon )%% Author :% Rafael Palac ios% Univers idad P o n t i f i c i a Comi l las% Madrid , Spain% Version : Apr /06 , Jun /06 , Aug /06 , Aug/06% Aug / 0 6 : f i x e d a problem ( found by Rodolphe Dewarrat ) r e l a t e d to southern% hemisphere coord ina tes .% Aug / 0 6 : cor rec ted m- L i n t warnings%- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -% Argument checking%e r r o r ( nargchk (2 , 2 , narg in ) ) ; %2 arguments requ i redn1= leng th ( Lat ) ;n2= leng th ( Lon ) ;i f ( n1~=n2 )e r r o r ( ' Lat and Lon vec to rs should have the same leng th ' ) ;end% Memory pre - a l l o c a t i o n%x=zeros ( n1 , 1 ) ;y=zeros ( n1 , 1 ) ;utmzone ( n1 , : ) = ' 60 X ' ;% Main Loop%209f o r i =1:n1l a =Lat ( i ) ;l o =Lon ( i ) ;sa = 6378137.000000 ; sb = 6356752.314245;%e = ( ( ( sa ˆ 2 ) - ( sb ˆ 2 ) ) ˆ 0.5 ) / sa ;e2 = ( ( ( sa ˆ 2 ) - ( sb ˆ 2 ) ) ˆ 0.5 ) / sb ;e2cuadrada = e2 ˆ 2 ;c = ( sa ˆ 2 ) / sb ;%alpha = ( sa - sb ) / sa ; %f%ablandamiento = 1 / alpha ; % 1/ fl a t = l a * ( p i / 180 ) ;lon = l o * ( p i / 180 ) ;Huso = f i x ( ( l o / 6 ) + 31) ;S = ( ( Huso * 6 ) - 183 ) ;de l taS = lon - ( S * ( p i / 180 ) ) ;i f ( la<-72) , Le t ra= 'C ' ;e l s e i f ( la<-64) , Le t ra= 'D ' ;e l s e i f ( la<-56) , Le t ra= 'E ' ;e l s e i f ( la<-48) , Le t ra= 'F ' ;e l s e i f ( la<-40) , Le t ra= 'G ' ;e l s e i f ( la<-32) , Le t ra= 'H ' ;e l s e i f ( la<-24) , Le t ra= ' J ' ;e l s e i f ( la<-16) , Le t ra= 'K ' ;e l s e i f ( la<-8) , Le t ra= ' L ' ;e l s e i f ( la<0) , Le t ra= 'M ' ;e l s e i f ( la<8) , Le t ra= 'N ' ;e l s e i f ( la<16) , Le t ra= 'P ' ;e l s e i f ( la<24) , Le t ra= 'Q ' ;e l s e i f ( la<32) , Le t ra= 'R ' ;e l s e i f ( la<40) , Le t ra= 'S ' ;e l s e i f ( la<48) , Le t ra= 'T ' ;e l s e i f ( la<56) , Le t ra= 'U ' ;e l s e i f ( la<64) , Le t ra= 'V ' ;e l s e i f ( la<72) , Le t ra= 'W' ;e lse Let ra= 'X ' ;enda = cos ( l a t ) * s in ( del taS ) ;eps i l on = 0.5 * log ( ( 1 + a ) / ( 1 - a ) ) ;nu = atan ( tan ( l a t ) / cos ( del taS ) ) - l a t ;v = ( c / ( ( 1 + ( e2cuadrada * ( cos ( l a t ) ) ˆ 2 ) ) ) ˆ 0.5 ) * 0.9996;ta = ( e2cuadrada / 2 ) * eps i l on ˆ 2 * ( cos ( l a t ) ) ˆ 2 ;a1 = s in ( 2 * l a t ) ;a2 = a1 * ( cos ( l a t ) ) ˆ 2 ;j 2 = l a t + ( a1 / 2 ) ;j 4 = ( ( 3 * j 2 ) + a2 ) / 4 ;j 6 = ( ( 5 * j 4 ) + ( a2 * ( cos ( l a t ) ) ˆ 2) ) / 3 ;a l f a = ( 3 / 4 ) * e2cuadrada ;beta = ( 5 / 3 ) * a l f a ˆ 2 ;gama = ( 35 / 27 ) * a l f a ˆ 3 ;Bm = 0.9996 * c * ( l a t - a l f a * j 2 + beta * j 4 - gama * j 6 ) ;210xx = eps i l on * v * ( 1 + ( ta / 3 ) ) + 500000;yy = nu * v * ( 1 + ta ) + Bm;i f ( yy<0)yy=9999999+yy ;endx ( i ) =xx ;y ( i ) =yy ;utmzone ( i , : ) = s p r i n t f ( '%02d %c ' ,Huso , Le t ra ) ;endD.25 setup matrix.md e f i n e g l o b a l ;ModelSetup ;LoadDataFi les ;number ce l ls = 3524;%leng th ( w e t t e d p a r i m e t e r s f r o m f i l e ) ;t o t a l l e n g t h = number cel ls ;Q = 100; %Cubic meters per secondus e leva t i on = 231;ds e leva t i on = 229;channe l mat r i x = zeros ( number cel ls , ma t r i x w id th ) ;%Def ine c e l l leng thschanne l mat r i x ( : , p length ) = t o t a l l e n g t h / number ce l ls ;%Def ine E leva t ionschanne l mat r i x ( : , pe leva t i on ) = bed e leva t ion ( 1 : number ce l ls ) ;%def ine depthschanne l mat r i x ( : , pdepth ) = NominalDepth ;s lope = ( us e leva t ion - ds e leva t i on ) / ( number ce l ls ) ;%channe l mat r i x ( : , pdepth ) = us e leva t i on - changes+6 - bed e leva t ion ( 1 :number ce l ls ) ;%def ine sec t ion c lasschanne l mat r i x ( : , pc lass ) =2;%def ine f r i c t i o n c o e f f i c i e n t s . Set to d e f a u l t value ;channe l mat r i x ( : , pk f ) = 0 .02 ;%def ine f low ra teschanne l mat r i x ( : , pQ) = Q;211f o r se tup index = 1: leng th ( channe l mat r i x )[A , Rh, Wp] = s e c t i o n p r o p e r t i e s ( channe l mat r i x ( setup index , pc lass ) ,channe l mat r i x ( setup index , pdepth ) , se tup index ) ;channe l mat r i x ( setup index , pArea ) = A;channe l mat r i x ( setup index , pRh) = Rh ;channe l mat r i x ( setup index ,pWp) = Wp;channe l mat r i x ( setup index , pX) = setup index ;%channe l mat r i x ( setup index , pdepth ) = us e leva t i on - setup index * slope + 4.4- channe l mat r i x ( setup index , pe leva t i on ) ;endload MeasuredDepths . matf o r i =1: leng th ( channe l mat r i x )channe l mat r i x ( i , pMeasDepth ) = i n t e r p 1 ( r ( : , 6 ) , r ( : , 3 ) , i ) ;endc l ea r setup index ;%def ine V e l o c i t i e schanne l mat r i x ( : , pve l ) = channe l mat r i x ( : , pQ) . / channe l mat r i x ( : , pArea ) ;%def ine To ta l Headschanne l mat r i x ( : , pTotalHead ) = channe l mat r i x ( : , pdepth ) + channe l mat r i x ( : ,pe leva t i on ) + channe l mat r i x ( : , pve l ) . ˆ 2 . / ( 2 . * g ) ;channe l mat r i x= reca lcu la te channe l params ( channe l mat r i x ) ;D.26 ProcessResultsMat.mPower = [ ] ;HeadChange = [ ] ;TurbArea = [ ] ;entranceDepth = [ ] ;f l ow = [ ] ;f o r i = 1 : s ize ( ResultsMat , 3 )Power = [ Power ; sum( ResultsMat ( : , pTurbA , i ) . * ResultsMat ( : , pvel , i ) . ˆ 3 * 0 . 5 * . 2 5 *rho ) ] ;HeadChange = [ HeadChange ; ResultsMat (1 , pdepth , i ) - ResultsMat (1 , pdepth , 2 ) ] ;TurbArea = [ TurbArea ; sum( ResultsMat ( : , pTurbA , i ) ) ] ;entranceDepth = [ entranceDepth ; ResultsMat (1 , pdepth , i ) ] ;f l ow = [ f low ; ResultsMat (2 ,pQ, i ) ] ;endD.27 depth sensitivity.m%Depth S e n s i t i v i t yModelSetup212NominalDepth = 4 . 0 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf2 = Performance ;CM2 = channe l mat r i x ;NominalDepth = 4 . 4 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf3 = Performance ;CM3 = channe l mat r i x ;NominalDepth = 5 . 0 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf8 = Performance ;CM8 = channe l mat r i x ;NominalDepth = 5 . 4 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;213channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf4 = Performance ;CM4 = channe l mat r i x ;NominalDepth = 6 . 0 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf5 = Performance ;CM5 = channe l mat r i x ;NominalDepth = 6 . 4 ;se tup mat r i xf p r i n t f ( ' So lv ing f o r k f ' ) ;k f = SolveKf ( NominalDepth ) ;channe l mat r i x ( : , pk f ) = k f ;f p r i n t f ( ' Running range of t u r b i n e constants ' ) ;r a n g e o f k t ;Perf6 = Performance ;CM6 = channe l mat r i x ;D.28 smooth.mf u n c t i o n out = smooth ( in , width )i f (mod( width , 2 ) ==0)width = width - 1 ;endi f ( s i ze ( in , 2 )>1)214f p r i n t f ( 'WARNING: Smooth should only be used wi th a 1D ar ray .\n ' ) ;endband = ( width +1) / 2 ;out = [ ] ;%Expanding band f o r the f i r s t measurementsf o r i = 1 : bandout = [ out ; mean( i n ( 1 : ( 2 * i - 1 ) ) ) ] ;endf o r i =(band+1) : ( leng th ( i n ) -band )%[ i , ( i - band ) : ( i +band ) ]out = [ out ; mean( i n ( ( i - band ) : ( i +band ) ) ) ] ;end%Col laps ing band f o r the l a s t measurementsf o r i = ( leng th ( i n ) -band+1) : leng th ( i n )out = [ out ; mean( i n ( ( i - ( 2 * ( leng th ( i n ) - i ) ) ) : l eng th ( i n ) ) ) ] ;end215

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

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

Comment

Related Items