A modeling study on tides in the Port of Vancouver

In the present paper, tides in the port of Vancouver Harbour have been investigated with a high-resolution three-dimensional hydrodynamic model based on FVCOM (Finite Volume Community Ocean Model). The model was evaluated against field observations including tidal elevations and tidal currents, and the evaluation showed that the model was in good agreement with the observational data. Using the model, we first investigated the horizontal distributions of tides, tidal currents, and tidally induced residual circulation, and then investigated the tidal asymmetry and dynamic mechanisms of tidal flows in the harbour. The tidal residual circulation shows a strong spatial pattern, which is associated with the local coastlines and variation of topography. The tidal asymmetry in the harbour is caused by different mechanisms, not only including the traditional factors, such as residual flows, the interaction between M2 and its overtide M4, but also the interaction of principal astronomical tides of O1, K1, and M2. The momentum balance is dominated by terms of the advection and the pressure gradient in First Narrows and Second Narrows, whereas terms of the local acceleration and the Coriolis are also important in the central harbour. The spatial variations of the momentum terms are strongly associated with the local changes in coastline and topography.


Introduction
The Port of Vancouver is located in the city of Vancouver and connected with the Pacific Ocean via Strait of Georgia ( Fig. 1) and is the largest and busiest harbour in Canada in terms of tonnage (McLaren 1994). The coastline inside the port has features with estuaries, headlands, islands, and infrastructure, such as terminals, docks, and bridge piers, which have spatial scales varying from several meters to ∼10 km. The distribution of the water depth in the harbour shows a strong spatial variation as well, from tidal flat areas to about 60 m at the central harbour. The harbour mainly consists of three parts, from west to east, First Narrows, the central harbour, and Second Narrows (Fig. 1). The central harbour shows an elliptic shape with a long axis of about 5 km in the east-west direction and a short axis of about 2.5 km in the north-south. First Narrows is a tidal channel connecting the central harbour with English Bay, while Second Narrows connects the central harbour with Indian Arm (Fig. 1). The mean width is about 880 m for First Narrows, and 330 m for Second Narrows. Tides are a major dynamic process in the Port of Vancouver. The tidal current enters the Strait of Georgia through the south passages, propagates into English Bay, and then enters into the harbour through First Narrows. The tides are mixed and dominated by the semi-diurnal M 2 component (12.42 h period) and the diurnal K 1 component (23.93 h period). Tidal currents in the two narrows are strong, with a magnitude up to 2 m s −1 (Isachsen and Pond 2000;Stacey and Pond 2003). The tidal process plays an important role in offshore infrastructure design, sediment transport, navigation safety, and biological processes in the harbour. So, understanding the hydrodynamic processes induced by tides in the harbour is important for the regular operation of the harbour, public safety, and marine environmental management.
Previous studies indicated that currents in the two narrows are dominated by tidal currents and the contributions of the wind forcing and freshwater are relatively weak Pond 1987, 1989). Using sea level and current data, de Young and Pond (1989) investigated the barotropic energy flux across the two narrows, and found that about 18% of the energy entering First Narrows was lost due to friction, mainly in the two narrows. Since the 1990s, a number of numerical models have been developed for the harbour for various purposes (e.g., de Young and Pond 1987;Stacey et al. 2002;Stacey andPond 2003, 2005). For instance, Stacey et al. (2002) developed a laterally averaged two-dimensional model, which covered the Burrard Inlet and Indian Arm with a high vertical resolution (1-3 m) and a horizontal resolution varying from 90 to 670 m. The model explained the propagation mechanism of internal tides in the Indian Arm and its influence on the exchange of the deep water. Li and Hondgins (2004) developed a three-dimensional (3-D) circulation model to study the effluent dispersion from a sewage outfall at First Narrows. Their model included tides, winds, and density stratification, and used a structured model mesh with 13 layers in the vertical direction and a uniform resolution of 75 m in the horizontal.
The complicated topography and coastlines in the harbour require a high-resolution model that is able to describe flow structures with small spatial scales in the harbour (Larson et al. 2003;Li and Hondgins 2004). Moreover, a model with enough resolution to represent the abrupt changes in the geometry and topography is necessary because the small-scale variations of the coastline and topography are known to lead to submesoscale motions characterized with spatial scales of 100 m -10 km (Hannah and Wright 1995;Zhao et al. 2006;Lynge et al. 2010;Poje et al. 2010;Xu and Xue 2011;Poje et al. 2014;McWilliams 2016). However, high-resolution models, with a spatial resolution to fit the spatial scales in the harbour, are still not available, though the computational power has been significantly increased in recent years. Recently, Fisheries and Oceans Canada plans to develop a national operational system based on the state-of-the-art hydrographic data and models with an aim to provide timely hydrographic information to vessels when they navigate in or approach ports. So, in this study, a model with an unstructured grid system is developed. The unstructured mesh allows us to use very high resolutions in areas of interest, and to fit complex coastline and topography. The arrangement of the paper is as follows: In Sect. 2, the model configuration is described. In Sect. 3, we present the model evaluation against the field observations, including water elevations and water currents. In Sect. 4, we first present the model results of tides, tidal currents, and the tidally induced residual circulation, and then we investigate the tidal asymmetry and momentum balance. Conclusions are given in Sect. 5.

Methods
The model used in this study is the Finite-Volume Community Ocean Model (FVCOM), which is a 3-D, finite-volume, unstructured grid ocean model (Chen et al. 2003(Chen et al. , 2007. The model is based on the triangle mesh system, which allows users to exploit variable resolution to fit complex coastline and topography. The model has a free surface, uses sigma coordinates, employs split-mode time stepping technique, and solves prognostic 3-D primitive equations. The model also contains an embedded second-order turbulence closure scheme, uses a Smagorinsky diffusivity to parameterize the horizontal eddy viscosity, and uses a wet-dry treatment to simulate flooding-drying process in tidal flats. The model domain covers a part of Strait of Georgia, Howe Sound, Burrard Inlet, Indian Arm, Pitt Lake, and Fraser River ( Fig. 2A). The coastline and bathymetric data were obtained from Canadian Hydrographic Service (CHS). The model uses 46806 nodes and 85804 triangular elements with size varying from ∼1 km in Strait of Georgia to ∼1 m in the harbour (Fig. 2B). In Second Narrows, for instance, the minimum element size is about 2 m to represent the effect of bridge bases. There are 21 sigma levels in the vertical direction with enhanced resolutions in the surface and the bottom layers. In the present paper, the model was run in a barotropic mode initialized using a constant temperature (18°C) and salinity (35) field. The open boundary forcing uses the amplitudes and phases of the tidal elevations derived from solutions of Foreman et al. (2000), a high-resolution, barotropic, tidal model with assimilating satellite altimeter data in the northeast Pacific Ocean. A sponge layer is utilized at the open boundary to reduce spurious barotropic resonant modes and internal model reflection within the model domain. The external time step is 0.2 s with a split number of 1.

Model evaluations
The model is evaluated against observational tidal elevation and tidal current data.

Tidal elevation evaluation
The hourly water elevation data used for the model evaluation are from two tidal gauges, 7735 (Vancouver) and 7795 (Point Atkinson), maintained by CHS. Gauge 7735 is located at the inner harbour, while 7795 is located at the northwest corner of Burrard Inlet (Fig. 1). Only data from February and March 2017 are used to keep consistency with the time of the observational current data. The modeled water elevations at the nodes nearest to locations of the tidal gauges are used. In the model-data comparison, two steps are used. We first compare the harmonic constants of specific tidal constituents and then compare the time series of tidal elevations to show the performance of the model on short timescales. The harmonic constants of the amplitudes and phases from the model and observations were computed using T-Tide (Pawlowicz et al. 2002), which is based on the theoretical analysis of Foreman (1977Foreman ( , 1978, and the results of six major tidal constituents are listed in Table 1. At the tidal gauge 7735, the tides are dominated by M 2 and K 1 with amplitudes of 0.94 and 0.86 m, respectively. The modeled results show similar amplitudes of 0.93 m for M 2 and 0.85 m for K 1 . The observed tidal form number (K 1 + O 1 )/(M 2 + S 2 ) is 1.12, which agrees well with the observed value of 1.13. The modeled phase of M 2 is 35.8°, slightly lower than that of observation by 0.2°. For K 1 , the modeled phase is lower than that of the observed result by 2.7°. For other tidal constituents, the differences between the model and observations are less than 0.01 m in amplitude and less than 3.7°in phase. For the tide gauge 7795, similar performance of the model was found. The model-data comparisons of the timedependent tidal elevations at the two tidal gauges are plotted in Fig. 3. It is worth noting that only the tidal component of the observed water elevations is used. To provide a quantitative comparison, we calculated four assessment parameters between the model results and observations following Zhong and Li (2006) and the parameters are the root mean square of the difference (rms), relative error (E), correlation coefficient (r), and skill (S). The definitions of the parameters can be found in Zhong and Li (2006). The skill (S) is the primary indicator of the goodness of fit between the model results and observations. Better fits are indicated by higher values, with a value of 0.5 indicating fair agreement. As we can see, the modeled tidal elevations follow the observations very closely at both stations. For the two  tidal gauges, the root mean square of the difference between the model and the observation is 0.05 m and the correlation coefficient and the skill are about 1.00.

Tidal current evaluation
Two kinds of current data are used in the model evaluation. One is the current data from four bins across Second Narrows collected with a horizontal Acoustic Doppler Current Profiler (h-ADCP) during the period of 2011-2012. The bins are located in the main channel of the Narrows with an interval of 20 m. The time interval of the velocity data is 15 min and the vertical location is 3 m below the water surface. The other kind is the current data collected in February-March 2017 by a ship-mounted ADCP for sections in First Narrows, the central harbour, and Second Narrows (locations are shown in Fig. 1).
As an example, we plotted the time series of the velocity data at bin 2 in Fig. 4. The currents are dominated by tides and show spring-neap cycles. The magnitude of the component u is much stronger than the component v. There are gaps in the data collection, especially in fall of 2011 and summer of 2012. In the model-data comparison, we use data in the period of February-March 2012, which has better data quality and moreover the stratification is weaker compared to other periods. The modeled depth-averaged currents at elements nearest to the bins are used. Tidal ellipse parameters from the model and observations of M 2 and K 1 are listed in Table 2. For M 2 , the modeled mean semi-major axis at the four sites is 1.47 m s −1 , which is slightly lower than the observed counterpart of 1.51 m s −1 ( Table 2). The semi-minor axis of both the modeled and observed M 2 tide rotates in the anti-clockwise direction with a mean difference of 0.001 m s −1 . The mean inclination from the model is 2.4°, close to the observed value of 2.5°. The mean difference of the M 2 phase between the model and the observations is 2.2°, corresponding to about 4.5 min. For K 1 , a similar model performance to that of M 2 can be found. Overall, the comparisons of the  (Table 2). At the same time, we can also see that the model slightly underestimates the tidal currents, especially the K 1 constituent.
The ship-mounted ADCP data were collected from 28 February to 4 March 2017 along 89 transits. For the present analysis, current data were vertically averaged into 50 m horizontal bins along the transits. Unlike the mooring data used above, the non-tidal component is not removed due to the short observation duration. Examples of the model and data comparison along selected transits are plotted in Fig. 5. The comparison shows that the modeled currents overall agree well with the observations, in both the flood and ebb directions. Both model results and observations show the strong currents in the two narrows, and in the central harbour a tidal dipole can be found with a cyclonic eddy on the north side of the main current and an anti-cyclonic eddy on the south side (Fig. 5C).
Comparisons of velocity components of all the ship-mounted ADCP data are plotted in Fig. 6. In First Narrows, the model agrees well with the observations in flood (u > 0 and v < 0). In ebb (u < 0 and v > 0) the model underestimates the velocity u (Fig. 6A), but overestimates the velocity v (Fig. 6B). The modeled current in Second Narrows matches well with the observed counterpart in ebb (u < 0), however, largely overestimates velocities in flood (u > 0, Figs. 6C and 6D). In the central harbour, the model is in good agreement with the observations in the u component, but in relatively poor agreement in the v component, due likely to the low magnitudes (Figs. 6E and 6F). To quantify the model performance, we analyzed the ship-mounted ADCP data using a quantitative method introduced by Hannah et al. (2001). In the method, four variables are calculated to quantify the modeled currents with the observed currents. The variables are a difference of the velocity magnitude (V d ), a difference of the velocity angle (Φ d ), the ratio of the kinetic energy in the velocity difference to the kinetic energy in the observed velocity (r) and the correlation coefficient (R) between the modelled and observed velocities. Among the four variables, r is the primary indicator of the goodness of fit between the model results and observations. Better fits are indicated by lower r values, with a value of 0.5 indicating a fair agreement. The definition of r is in which v m and v o are the modeled and the observed velocity vectors, respectively. The definition of the other three variables can be found in Hannah et al. (2001).
The comparison results are listed in Table 3. In First Narrows, modeled mean speed is 0.71 m s −1 , higher than that of the observed by 0.11 m s −1 . The model overestimates the mean speed in the central harbour with an increment of 0.03 m s −1 , but slightly underestimates the speed in Second Narrows. The differences in the velocity magnitudes in the three areas are 0.18, 0.12, and 0.26 m s −1 , and differences in the velocity angles are 26.3°, 47.9°, and 37.0°f or the three areas, respectively. The correlation coefficients are high in the two narrows; 0.97 for First Narrows, and 0.95 for Second Narrows, while they are lower in the central harbour at 0.70. The value of r is 0.09 for First Narrows, 0.08 for Second Narrows, and 0.53 for the central harbour. Overall, the model-data comparisons show that the model results are in good agreement with the in situ observations in First Narrows and Second Narrows, and fair agreement is found in the central harbour. For the two Narrows, the good agreement is expected to benefit from the strong tidal currents, which are much stronger than the non-tidal component caused by non-tidal forcing, such as winds and density gradients. For the inner harbour, however, the worse agreement is probably due to the simplification in the model and uncertainties of the model input. For example, effects of winds and the density gradient are not included in the present version of the model. At the same time, we recognized that the time duration of the ship-mounted ADCP data is too short to remove the non-tidal components, which are not included in the model.

Results
We will present the model results of the tidal elevations and tidal currents, the tidally induced residual circulation, tidal asymmetry, energy flux, and the momentum balance. Table 3. Goodness-of-fit statistics from the comparison between observed and modelled current velocities along the section shown in Fig. 1.

Tidal amplitudes and phases
The tidal amplitude of M 2 falls in a narrow range of 0.95-0.99 m in the harbour (Fig. 7A). Nevertheless, there are important variations, such as the amplitude decrease of about 2 cm at the western inlet of First Narrows. The magnitude then increases in the narrows with a spatial variation associated with the local coastline and topography, for instance, at the east outlet of First Narrows, where the tidal magnitude varies around the headland of Brockton Point. In the central harbour, the M 2 amplitude increases from west to east. Changes in M 2 phase mainly occur at the front of the two Narrows (Fig. 7B). The phase is about 35°at the west entrance of First Narrows, and it increases quickly to 45°in the end of the Narrows. The difference reaches 10°, corresponding to 20.7 min. A similar difference occurs from 45°to 55°in the east of the Narrows. This spatial change in phase is consistent with the observations of de Young and Pond (1989) and Stacey (2005). The amplitude of K 1 constituent ranges from 0.90 to 0.94 m in the harbour. Similar to M 2 , a sharp decrease of the K 1 amplitude can be found at the western inlet of First Narrows (Fig. 7C). A 2 cm decrease can also be found at Second Narrows. Similar to M 2 , the phase of K 1 increases from west to east in the harbour (Fig. 7D), a phase lag is about 10°in the two narrows, corresponding to 40 min.

Tidal currents
To show the tidal currents, we plotted the averaged tidal currents in flood durations and ebb durations instead of the tidal harmonic constants used for the tidal elevations. The duration of flood-ebb is defined as the time with the east-west current at the west inlet of First Narrows after a sensitivity test of the duration to different locations, for example, at the central harbour and the inlet of Second Narrows. During the flood duration, tidal currents flow from west to east. The speed of the tidal flow is usually <0.2 m s −1 in English Bay, but increases to 1.5 m s −1 in First Narrows (Fig. 8A). We can see clearly that the main current flows southeast along the channel of the narrows. At the same time, recirculation structures can be seen as well in the nearshore water along the narrows. Along the south coast, for example, a recirculation flows with an opposite direction to the main flow. Along the northern coast, a recirculation structure can be found with a lower magnitude and smaller spatial scale than the southern one. The formation of the recirculation structures is likely related to the changes in coastline, for example, the presence of Brockton Point, which forms a headland and separates the flood currents at the front of it, and causes the eddy along the south coast. In the central harbour, the current becomes narrower and weaker (Fig. 8A), extending to the southern coast, where it is separated into two branches. The main branch flows northeast entering into Second Narrows and part of it forms the cyclonic eddy covering most of the central harbour. The minor branch flows west and forms a small anticyclonic eddy covering the southwest water. Similar to First Narrows, a significant increase of the tidal flow speed can be found in Second Narrows with the speed of 2.5 m s −1 . The tidal current then decreases in the Indian Arm. Similarly to First Narrows, recirculation structures can also be found in Second Narrows in areas, where coastline and the topography change. In ebb durations, the current in Second Narrows generally flows from east to west with a magnitude of 2.5 m s −1 and then significantly decreases in the central harbour (Fig. 8B), where the ebb current mainly flows west and the main part of the current is along the northern part of the harbour. The water converges and enters into First Narrows, where the flow magnitude reaches 1.5 m s −1 . Similar to the flood, recirculation structures can be found in areas where abrupt changes in coastline and topography occur (Fig. 8B). The current patterns in flood and ebb durations are similar to the tidal current charts of Canadian Hydrographic Services (1973). The current moves from west to east in flood, while it moves from east to west in ebb. In the two narrows, the current is strong and the direction is overall along the narrows. In the central harbour, however, the current shows complex structures with the variations of the coastline and topography.
Our model results also show vortex-pairs at the eastern end of First Narrows after slack ebb tides. Based on previous studies, vortex-pairs, or vortex dipoles, are commonly found in strait-basin systems, usually forming at the strait outlet after slack tides due to the abrupt widening of the basin (Imasato 1983;Kashiwai 1984;Signell and Geyer 1991). The vortex pairs self-propagate away from straits into basins during each tidal cycle and have been considered to play an important role in exchanges of mass and momentum between basins and the open water (Fujiwara et al. 1994;Wells and van Heijst 2003;Nicolau del Roure et al. 2009). As an example, we plotted the vortex-pairs in a spring tidal cycle in Fig. 9. The intensity of the vortex is represented by the relative vorticity, defined as ζ = ( ∂v/ ∂x) -( ∂u/ ∂y), where u (v) represents the velocity component in the eastward (northward) direction. At 1.5 h after slack water of the ebb tide, the tidal velocity at the outlet of First Narrows quickly increases from zero to 0.8 m s −1 (Fig. 9A). One can see a cyclonic eddy at the north side of the outlet with a magnitude of vorticity of 4 × 10 −3 s −1 and an anti-cyclonic eddy on the south side of the outlet with a comparable vorticity magnitude. With the development of the flooding tide, the intensity of the dipole strengthens and the moving speed enhances until it meets the coast, where the dipole breaks (Figs. 9B-9D). The dipole breaks up at the front of the south coast. The south eddy was reflected to the southwest while the north eddy moves north and dissipates in the central harbour. The tidal dipole lasts 4-5 h. The dynamics of tidal dipoles have been extensively studied using both experimental and numerical methods (Imasato 1983;Kashiwai 1984;Signell and Geyer 1991;Wells and van Heijst 2003). The vortex-pair originates from the low pressure areas that are generated downstream behind a headland by the nonlinearity of the centrifugal force of the tidal current flowing with a large curvature through a narrow channel (Imasato 1983). The propagation of dipoles is usually related to parameters of the tidal inlet width, the maximum velocity through the inlet, the tidal period, and the channel length (Kashiwai 1984;Wells and van Heijst 2003;Bryant et al. 2012).

Residual circulation
Residual circulation caused by tides plays an important role in long-term transport of water mass, heat, salt, nutrients, and sediment (Nihoul and Ronday 1975;Tee 1977;Robinson 1981). In this study, the residual current is defined as the mean of the depthaveraged currents subtracting the tidal component, and the results are plotted in Fig. 10. In the western inlet of First Narrows, the residual circulation structure features an eastern current along the southern slope and a western current along the northern slope (Fig. 10). The eastern current is caused by the stronger flood currents, whereas the western current is caused by the stronger ebb currents according to the current patterns in flood and ebb ( Fig. 8). In First Narrows, the residual circulation is dominated by an anti-cyclonic eddy that covers the narrows. The eddy in the southern coast of First Narrows is likely associated with the presence of Brockton Point, which reflects the southern flow to the north at the front of the headland in fold phases and forms an anti-cyclonic eddy along the south coast of the Narrows (Fig. 8A). In the central harbour, the residual circulation is dominated by a cyclonic eddy, which covers most of the central harbour. The formation mechanism of the eddy can be related to the geometry of the harbour and variations of the coastline and topography, which lead to the different pathways of tidal flows in flood and ebb, and causes nonlinear interactions between flows and topography (Huthnance 1973;Loder 1980;Marinone 1997). Moreover, the tidal dipole from First Narrows is also responsible for the formation of the residual eddy according to previous studies (e.g., Tee 1977;Imasato 1983;Fujiwara et al. 1994). The residual circulation in Second Narrows also features eddies, which are strongly related to variations of coastline and topography (Fig. 10B).

Tidal asymmetry
Tidal current asymmetry, defined as the periodic differences in tidal velocities between flood and ebb duration, is an important metric for understanding the spatial variability in the underlying tidal dynamics, the net sediment transport, and morphological changes (Dronkers 1986;. Physically, the occurrence of the tidal current asymmetry is associated with irregular bathymetry, and interactions among tidal constituents (Wu et al. 2011;Jewell et al. 2012;Guo et al. 2016). A commonly known interaction is that a principal tide interacts with either itself or with its higher harmonic tides (over tides). For example, in semi-diurnal tidal regimes, tidal asymmetry can be produced by M 2 interacting with itself to produce a mean flow or by the phase relationship between M 2 and its first harmonic, M 4 , which is generated nonlinearly by the M 2 tide in shallow waters Speer and Aubrey 1985). In the absence of nonlinear dynamics, tidal current asymmetry can be produced through the interaction of the principal tides (Ranasinghe and Pattiaratchi 2000;Maren and Gerritsen 2012). For instance, in mixed, mainly semi-diurnal tidal regimes, tidal asymmetry is formed due to the combination of the lunar (K 1 ), lunisolar (O 1 ) diurnal tides with the semi-diurnal constituent M 2 (Woodworth et al. 2005;Nidzieko 2010;Jewell et al. 2012;Maren and Gerritsen 2012). Following Aldridge (1997), tidal asymmetry index in the present study is estimated by comparing the maximum flood and ebb velocities in tidal cycles: γ = max juj flood − max juj ebb max juj flood + max juj ebb (2) where max |u| flood and max |u| ebb are the maximum speed during the flood and ebb durations, respectively. For the tidal asymmetry index γ > 0, the tidal current is flood-dominant, while it is ebb-dominant for γ < 0. In this study, the flood (ebb) duration is defined as the time period when the mean velocity flows east (west) across a section at the western entrance to First Narrows (selected after considering sections in both First and Second Narrows). The tidal asymmetry shows a strong spatial variation (Fig. 11A), especially at the ends of the two narrows. For example, the negative index can be found in the west area of First Narrows because of the strong tidal jets during the ebb duration. On the contrary, at the area of the eastern end of the Narrows, the asymmetry index is positive with a magnitude of 0.6. A similar feature can also be found at the ends of Second Narrows. Besides the spatial variations around the two Narrows, strong tidal asymmetry can be found in the central harbour. Positive index occurs in the southern part of the central harbour (Fig. 11A), whereas negative index can be found in the northern part of the harbour.
According to previous studies, the tidal current asymmetry is usually caused by three processes, which include (i) interactions of M 2 and its first higher harmonic M 4 (M 2 -M 4 ), (ii) interactions of principal tidal constituents of O 1 , K 1 , and M 2 (O 1 -K 1 -M 2 ), and (iii) the residual flow (Woodworth et al. 2005;Nidzieko 2010). Using the modeled tidal flows, we next analyze the specific contribution of each process to the tidal asymmetry. To do so, the instantaneous velocity, u, is assumed to consist of residual flow and four tidal harmonic components: uðx,y,tÞ = u res ðx,y,tÞ + u m2 ðx,y,tÞ + u m4 ðx,y,tÞ + u o1 ðx,y,tÞ + u k1 ðx,y,tÞ (3) in which the time series of velocities at each model element is derived from the harmonic analysis using T-Tide. It is worthwhile to note here that the residual component is the results with the removal of main principal tides (Q 1 , O 1 , K 1 , P 1 , N 2 , S 2 , and M 2 ) and main higher harmonic tides (M 4 , M 6 , K 2 , and MK 3 ). The velocity of one component, for example, M 2 -M 4 , can be estimated by removing the other components (residual, O 1 and K 1 ) in eq. (3). The tidal asymmetry of the mechanism can be estimated by replacing the numerator of eq. (2). Because we aim to find specific contributions to the tidal asymmetry of the default case shown in Fig. 11A, the denominator of eq. (2) in the calculation of the specific contribution of one mechanism still use the values of the total currents (default case). To show the horizontal structures clearly, tidal asymmetry index was multiplied by 4 (Figs. 11B-11D). The contributions of the three mechanisms to the tidal asymmetry vary across the harbour. Residual flows dominate the asymmetry, particularly in the nearshore areas, where the coastline is complex and the topography is irregular (Fig. 11B). For example, in the northwest and southwest parts of the central harbour we can see clearly the effect of residual flow overweighs M 2 -M 4 and O 1 -K 1 -M 2 factors (Fig. 11B). In contrast, the interaction of M 2 and M 4 dominates the tidal asymmetry in the two ends of both narrows (Fig. 11C). The presence of M 2 -M 4 leads to ebb Fig. 11. Tidal asymmetry index due to all tides (A), residual (B), M 2 -M 4 (C), and O 1 -K 1 -M 2 (D). To show the horizontal structures clearly, tidal asymmetry indices were multiplied by 4 (Figs. 11B-11D). The stars indicate the three locations used in Fig. 12. dominance in the western ends and flood dominance in the eastern ends. The contribution of O 1 -K 1 -M 2 to the asymmetry mainly occurs in the southern part of the central harbour covering from the east off the outlet of First Narrows to the southeast coast (Fig. 11D).
To examine the different mechanisms that lead to the tidal current asymmetry, we plotted time series of the velocity u component at three selected sites (Fig. 12). The flood-or ebb-dominance by M 2 -M 4 interaction is considered to be controlled by the phase relationship between the M 2 and M 4 constituents given by θ = 2φ u (M 2 ) -φ u (M 4 ), where φ u (M 2 ) and φ u (M 4 ) are the phase of M 2 and M 4 currents, respectively (Ranasinghe and Pattiaratchi 2000;Maren and Gerritsen 2012). For a phase difference θ with cos θ > 0, the tide is flood-dominant asymmetric, whereas, for a phase difference θ with cos θ < 0, the tide is ebb-dominant asymmetric. The magnitude of the tidal asymmetry is considered to be proportional to the amplitudes of M 2 and M 4 . Similarly to M 2 -M 4 , the qualitative contribution of O 1 -K 1 -M 2 to tidal asymmetry can also be identified by the phase difference between them, which is defined as θ = φ u (O 1 ) + φ u (K 1 ) -φ u (M 2 ), where φ u (O 1 ), φ u (K 1 ), and φ u (M 2 ) are respective phases of O 1 , K 1 , and M 2 currents (Maren and Gerritsen 2012). For a phase difference θ with cos θ > 0, the tide is flood-dominant asymmetric, whereas, for a phase difference θ with cos θ < 0, the tide is ebb-dominant asymmetric. The magnitude of the tidal asymmetry is considered to be proportional to the amplitudes of O 1 , K 1 , and M 2 .
The time series of the residual current at site 1 is negative (Fig. 12A) with a mean magnitude of 0.19 m s −1 , indicating that the contribution of the residual is the ebb-dominant asymmetry. The tidal current due to M 2 and M 4 shows the ebb-dominant asymmetry as well. The phase difference of M 2 and M 4 is θ = 165.1°, which indicates that M 4 is almost in phase with M 2 in ebb, but out of phase in flood. The combination of the semidiurnal and diurnal tides produces a negative asymmetry index (Fig. 12G) with a phase difference of θ = 133.7°, but with a much weaker magnitude than those caused by M 2 -M 4 and residual currents. On the contrary to site 1, the time series at site 2 shows a flood-dominant asymmetry (Figs. 12B and 12E). The magnitude of the residual current is 0.09 m s −1 , about half of the difference between the magnitude of the flood and ebb of M 2 -M 4 , which has a phase difference of θ = 358.4°, indicating that M 2 and M 4 are in phase in flood, while out of phase in ebb. The diurnal inequality in O 1 -K 1 -M 2 is amplified more in the flood durations with the phase difference of θ = 49.5°, however, the intensity of asymmetry is much weaker (Fig. 12H). Site 3 shows a flood-dominant asymmetry primarily controlled by the residual and O 1 -K 1 -M 2 , rather than the M 2 -M 4 interaction. The magnitude of residual and O 1 -K 1 -M 2 reaches 0.25 m s −1 , which is more than five times of that of M 2 -M 4 despite an M 2 -M 4 phase difference of θ = 22.1°.
Overall, the analysis above indicates that the tidal asymmetry in the harbour is flood-dominated and caused by different mechanisms, which not only include the traditional factors, the residual flows, and the interaction between M 2 and its overtide M 4 , but also the interaction of principal astronomical tides of O 1 , K 1 , and M 2 . The contribution of each factor shows a strong spatial variation, which is caused by abrupt changes in shorelines and topography. The influence of the tidal asymmetry on the sediment transport is important, however, the detailed discussion about the influence mechanisms will be presented in a separate paper.

Energy flux
As the barotropic tides propagate over tough topographic features, a portion of the barotropic energy is lost due to local dissipation by the bottom friction, horizontal diffusion, and conversion to internal tides (Seim et al. 2006;Easton et al. 2012;Kang and Fringer 2012). In Vancouver Harbour, the energy loss of barotropic tides mainly occurred in First and Second Narrows based on an investigation of the observational water elevations at three points Pond 1987, 1989;Stacey 2005). The two Narrows dissipated about 18% of the incoming barotropic energy and each of them was responsible for the same amount of 8%. De Young and Pond (1989) reported that the bottom friction dominated the energy dissipation mechanisms. In this study, we use the model results to show a more detailed distribution of the energy loss in First and Second Narrows. Following Zhong and Li (2006), tidal energy loss in a controlled volume can be calculated by the total net tidal energy across the boundary of the controlled volume, and the net tidal energy along the unit width is written as where η is the water elevation, u is the depth-averaged horizontal velocity vector, h is the total water depth, ρ 0 is the water density, and g is the acceleration of gravity. The two terms on the right-hand side of eq. (4) represent the energy flux due to potential energy and kinetic energy, respectively. We first interpolate the water elevations at model element nodes into the element centers and then average the sum of the two terms over the model run duration of 28 days. Using the model results, the total energy flux is shown in Fig. 13. The energy flux shows a strong spatial variation in magnitudes and directions. A large incoming tidal energy flux is seen to flow across the west inlet of First Narrows. The flux flows along the Narrows from west to east in the main Narrows. In the southern side of the narrows, the net energy flux can be seen to flow from east to west, however, with a much weaker magnitude compared to that in the main Narrows. A similar spatial feature can be seen in the Second Narrows, where the energy flux moves from the west along the main Narrows with local structures caused by the variations in topography and coastlines. Based on the definition of eq. (4), the total energy loss in the Vancouver Harbour can be calculated by integrating the net tidal flux across the section normal to flow at the west head of First Narrows. The energy loss over the area east to the section is 12.3 MW. Similarly, net energy flux based on eq. (4) of tidal constituents of M 2 and K 1 are also calculated. The horizontal features of the energy flux of the two constituents are very similar to Fig. 13, but the magnitude of K 1 flux is much smaller (6-7 times) than that of M 2 flux (not shown). The net energy through the west head of First Narrows is 8.7 MW for M 2 and 1.3 MW for K 1 , which are quite close to the results of de Young and Pond (1987), who reported that there is 7.1-8.1 MW lost from M 2 and 1.0 MW from K 1 .

Momentum balance
To understand the dynamics of the tidal flows, in this section, we examine the dominant terms of the momentum equations. Following Hench and Luettich (2003), the terms in the momentum equations were rotated from x-y components into a local streamwise-normal (s-n) coordinate system. In the s direction, the momentum is balanced between the terms of local streamwise acceleration, streamwise accelerations, streamwise pressure gradient, and bottom friction. In the n direction, the terms in the balance include the local rotary acceleration, centrifugal acceleration, Coriolis acceleration, and normal direction pressure gradient force. Following Hench and Luettich (2003) the final form of the equations in the s-n coordinate system are  where U s is the streamwise velocity, α is the streamline angle (the angle between the positive x-axis and the local flow vector), R s is the streamwise radius of curvature, and τ b is the bottom shear stress. H is the total water depth.
Because the model output is based on the triangle mesh system, we first linearly interpolated the time series of model results into a regular mesh system to solve the derivative terms conveniently. The standard deviations (SD) of each term in eq. (5) based on the model run in a period of 28 days are plotted in Fig. 14. Momentum terms in the streamwise direction show a strong spatial variation, and the magnitudes of the terms in the two Narrows are much higher than those in the central harbour (Figs. 14A-14D). Compared to First Narrows, the magnitudes in Second Narrows are relatively higher. In the two Narrows, the momentum balance in the streamwise direction is dominated by the advection term (Fig. 14B) and the pressure gradient term (Fig. 14C), and they are higher than those of the local acceleration (Fig. 14A) and bottom friction terms (Fig. 14D) by 1-2 orders of magnitude. This balance feature indicates that the magnitudes of the advection term and the pressure gradient are comparable, but signs are opposite. In the two Narrows, the spatial distribution of the advection and the pressure gradient is associated with the variation of the coastlines and topography. The predominant magnitude of the streamwise acceleration and the pressure gradient occurs around domes and headlands (Figs. 14B and 14C), where the higher pressure gradient speeds up the flow velocity in the downstream direction. This dynamic mechanism is commonly found in tidal straits with variable coastline and topography (Hench and Luettich 2003). The magnitude of the term of the local acceleration is comparable with that of the bottom friction term. In contrast, in the central harbour, the momentum is mainly balanced among the local acceleration, advection, and pressure gradient. The bottom friction in the inner harbour is lower mainly because of the deeper water depth and the relatively lower tidal velocity.
In the rotary (or normal) direction, the momentum balance is dominated by the centrifugal acceleration (Fig. 14F) and the pressure gradient (Fig. 14G) almost everywhere in the harbour, especially in the two Narrows, where these two terms are larger than the rotary acceleration and the Coriolis acceleration by 1-2 orders of magnitude (Figs. 14F and 14H). This means that the rotation of the flow (or vorticity) is also dominated by the advection and the pressure gradients, which is strongly associated with the variations of local coastlines and topography. The rotary acceleration and the Coriolis acceleration are comparable in magnitude. The Coriolis term nicely maps out the locations of the tidal jets in the inner harbour where the tidal currents are relatively stronger (Fig. 14E).
To have a quantitative comparison, we list the standard deviation of each term of eq. (5) in Table 4. At the same time, correlation coefficient between the pressure gradient and other terms are also calculated. In First Narrows, the standard deviation of the streamwise advection term is about 70% of the pressure gradient, the local acceleration is about 33%, and the bottom friction is about 35%. The correlation coefficient between the pressure gradient and the advection is 0.90 and bottom friction is 0.63. In the normal direction, the standard deviation of the advection term is about 71% of the pressure gradient, and the Coriolis force is 70% as well. This means that both the Coriolis term and the rotary advection are equally important for the vorticity formation in First Narrows. The local rotary acceleration is relatively small, about 23% of the pressure gradient. The correlation coefficient of the advection and Coriolis are 0.76 and 0.68 of the pressure gradient term, higher than the local rotary acceleration of 0.15. In the central harbour, the standard deviation of the streamwise local acceleration is about 80% of the pressure gradient, which is 1.5 times the bottom friction, and 3 times the streamwise advection. The correlation coefficient for the local acceleration is 0.83, the bottom friction is 0.54, and the advection term is 0.25. In the normal direction, the pressure gradient is mainly balanced by the advection term and the Coriolis term, which are equally important in the central harbour. The standard deviation of the advection and Coriolis force is 50% of the pressure gradient, and the local rotary acceleration is 2%. The correlation coefficient of Coriolis reaches 0.90, and the advection reaches 0.84. Similar to First Narrows, the pressure gradient in the streamwise direction in Second Narrows is by the horizontal advection, which is about 59% of the pressure gradient, and the bottom friction and local acceleration is 35% and 27%, respectively. The variation of the bottom friction and the advection term are highly correlated to the pressure gradient, with coefficients of 0.8 and 0.94, respectively. In the rotary direction, the standard deviation of the advection term is about 69% of the pressure gradient, Coriolis is 34%, and local rotary acceleration is 5%. Both the advection term and Coriolis force are highly correlated to the pressure with a coefficient of 0.84. Overall, the momentum balance is mainly dominated by the advection and the pressure gradient terms in the two Narrows, whereas the mainly balanced terms of the local acceleration and the Coriolis force are also important in the central harbour.

Conclusions
In the present paper, the tidal processes in Vancouver Harbour have been investigated with a high-resolution 3-D hydrodynamic model based on Finite Volume Community Ocean Model (FVCOM). The model uses the flexible triangle mesh system with a varying model resolution from a couple of kilometres in the offshore waters to several metres in the near-shore waters. The model was evaluated against observed tidal elevations and tidal currents and the evaluation showed that the model was in good agreement with the observations. The model shows that tides in Vancouver Harbour have a strong spatial variation, which is mainly dominated by the two strait-bay systems associated with First Narrows and Second Narrows. The tidal amplitude of M 2 shows a small decrease at the entrance of First Narrows, but then overall increases in the central harbour and Second Narrows. In contrast, the K 1 consistently decreases from First Narrows to Second Narrows. The phase of tides increases from the west to east, and significant decreases mainly occur at the entrances of the two Narrows, with about 10°for both M 2 and K 1 . Similarly to the tidal elevation, tidal current is controlled by the two strait-bay systems as well. The tidal speed is high in the two Narrows, and the current fields are strongly affected by the abrupt coastlines and irregular topography. Using the model output, we further analyzed tidal asymmetry and the tidal dynamics by comparing the momentum balance. The analysis indicated a strong tidal asymmetry in the harbour caused by not only the traditional factors, the residual flows, and the interaction between M 2 and its overtide M 4 , but also interaction of principal astronomical tides of O 1 , K 1 , and M 2 . The magnitudes of the momentum terms in the two Narrows are usually higher than those in the inner harbour. The momentum balance is mainly dominated by the term of the advection and the pressure gradient in the two narrows, whereas the terms of the local acceleration and the Coriolis force are also important in the central harbour. Overall, in the present study a high-resolution model was developed. Compared to the previous studies, the model has enough resolution to represent the abrupt changes in the geometry and topography of the harbour and to simulate submesoscale motions in the harbour. The model can be used to investigate the physical, environmental issues, and biological processes, for example, the effect of the tidal dipole on the heat transport and the sediment erosion and deposition in the ambient water.