Circulations in the Pearl River Estuary: Observation and Modeling

Abstract

This chapter reports a cruise survey on the Pearl River Estuary (PRE) and adjacent costal water in the period between May 3, 2014 and May 11, 2014. The circulation and salinity structure were sampled for different tidal phases. With the cruise data, a “sandwich” structure of the lateral salinity distribution and a two-layer structure of longitudinal circulation were identified, together with high variations influenced by wind and tide. Furthermore, longitudinally orientated convergence or divergence of the lateral velocity close to the channel location for certain tidal conditions was observed. The finite volume community ocean model (FVCOM) is configured and run with high spatial resolution of 100 m in the PRE. An atmospheric model, the Weather Research and Forecasting (WRF) Model, is also run to provide high spatial and temporal resolution of atmospheric forcing for the FVCOM. The FVCOM modeling skill assessment is conducted using the cruise salinity and velocity data, as well as water levels, showing that the model can well simulate the velocity and salinity structures. The numerical model reveals that there is a strong neap-spring cycle for the PRE de-tided circulation with 0.37 m s−1 during the neap tide about 42% stronger than that (0.26 m s−1) during the spring tide in the surface layer.

1. Introduction

The Pearl River Estuary (PRE) is located along the coast of Guangdong Province in China between 113°E and 115°E and 21°N and 23°N, in a trumpet-like shape with a width of ∼5 km at the northern (upper) end and ∼35 km at the southern (lower) end ( Figure 1 ). The Pearl River freshwater flows into the northern South China Sea (NSCS) through eight inlets, named Modaomen, Humen, Hongqili, Jiaomen, Jitimen, Hengmen, Yamen, and Hutiaomen, with an annual average discharge of 10,000 m³ s⁻¹ [1]. The river discharge reaches a maximum of 20,000 m³ s⁻¹ in the wet summer and has a minimum of 3600 m³ s⁻¹ in the dry winter, combining the discharge of all eight inlets together. In particular, four major inlets (Humen, Jiaomen, Hongqimen, and Hengmen) located along the west shore of the PRE contribute more than half of the total discharge [2]. There are two longitudinal channels in the central region and on the eastern side, connecting the Pearl River Estuary with the NSCS. The east channel connects the coastal water through the Urmston Road near the east bank of the PRE and the west one through the Lantau Channel along the middle estuary. Except for these channels where the water depth varies between 5 and 15 m, most of the PRE is shallow with a water depth between 2 and 10 m; outside the estuary, the water depth increases to more than 20 m, and the isobaths are approximately parallel with the coastline with a strong cross-shelf gradient on the shelf. Furthermore, the coastline is complex around the PRE and its adjacent shelf, with many small islands scattered around the estuary mouth and coastal bays. These geometric and topographic features can greatly affect the circulation in both the estuary and the coastal region.

Figure 1.

The northern South China Sea region (a) and the Pearl River estuary (b) with bathymetry contours (meter).

Generally, the circulations of the PRE consist of the buoyant gravitational circulation associated with the freshwater discharge and the complicated bathymetry, the tidal currents, and the wind-driven currents [1, 3]. The PRE is a micro-tidal estuary, and the M2 and K1 constituents dominate the tidal currents, which have a 1.0 m tidal range inside the PRE [4]. The mean tidal range just outside the PRE is 0.85–0.95 m [2]. The tidal currents also can change stratification structures in the PRE and play a vital role in the vertical mixing, especially in the shallow near-shore region [5, 6]. Using the Princeton Ocean Model (POM) with curvilinear orthogonal grid to study the Pearl River plume, Wong et al. [3] found that the water column near the head of the estuary becomes well mixing owing to tidal mixing and the moderate southwesterly wind has limited effects on the circulation inside the PRE [3, 7]. The circulation in the upper estuary and northern part of the lower PRE is characterized by the seaward gravitational current. However, the southern part of the lower estuary is governed by the competition between gravitational current and geostrophic intrusive currents from the shelf [5]. Based on cruise survey data and numerical modeling of the PRE for the summer time, Pan et al. [6] revealed the role that wind and tidal forcing play in determining plume dynamical properties and concluded that the density structure in the PRE can be modulated by the wind with strong stratification under the southwesterly wind and weak stratification under the southeasterly wind [8]. Under the southeasterly wind during spring tide, an anticyclonic circulation bulge appeared, with a Kelvin number of 0.93 and a supercritical Froude number of 1.45. Under southwesterly wind during a neap tide, the plume displayed a diffusive front, with a Kelvin number of 2.9 and a Froude number of only 0.07 [8]. The turbulent mixing in PRE is governed by the tidal strength, and the wind stirring tends to increase the internal shear instability mixing in the recirculating plume bulge [6]. Observation and modeling reveal that strong mixing appeared in the bottom layer on the larger ebb, caused by the vigorous bottom stress.

The spreading of the plume in the surrounding coastal waters is controlled by the wind stress and ambient currents. The easterly or southeasterly wind drives the plume westward, and southerly or southwesterly wind forces the plume eastward [3, 8, 9]. The wind-driven coastal current closely interacts with the estuarine circulation in the lower part of the PRE [1]. Ou used the field measurement data to investigate the dynamics of the buoyant plume near the PRE, suggesting four major horizontal plume types in summer [9]. These types are offshore bulge spreading, westward alongshore spreading, eastward offshore spreading, and symmetrical alongshore spreading. In this chapter, we report detailed cruise observations on the PRE water current and the salinity structure. A numerical model of high spatial resolution is developed in the PRE, which can well simulate the circulations of the PRE.

2. Cruise observation

A hydrographic cruise survey for collecting in situ data of marine environmental parameters in the PRE was conducted during the period between May 3, 2014 and May 11, 2014. The instruments used to collect the in situ data in this cruise survey were a sea surface sampling Seabird CTD (conductivity, temperature, and depth) for sea surface sampling, an RDI 600 kHz acoustic Doppler current profiler (ADCP), and a Sea and Sun Technology CTD. The RDI ADCP was mounted on the right side of the cruise ship. With these instruments, the underway surface salinity, temperature, density, and underway profiles of currents were collected, and additionally, the CTD profiles at a number of CTD cast stations along the cruise transects were also collected. Figure 2a displays the tidal conditions at the Shek Pik Station, located in the south of the Lantau Island near the mouth of the PRE (as shown in Figure 3c ), and the observed periods of different transects. The 6-hourly wind conditions from the NCEP/NCAR reanalysis data over the PRE region during the cruise period are shown in Figure 2b . The cruise transects and CTD cast stations across the PRE are illustrated in Figure 3 . During this cruise, a typical East Asian monsoon transition process happened with relatively weak wind strength. The northeasterly/easterly winds prevailed during a spring to neap tide transition from May 4, 2014 to May 6, 2014 ( Figure 3a ). From May 7, 2014 to May 9, 2014, the winds turned to be the easterly/southeasterly and the spring tide changed to neap tide ( Figure 3b ), and the southerly/southwesterly winds prevailed during a neap tide from May 10, 2014 to May 11, 2014 ( Figure 3c ).

Figure 2.

Time periods of cruise transects together with the sea level at the Shek Pik Station (the location is shown in Figure 3c ) (a). The 6-hourly wind condition from the NCEP/NCAR reanalysis data over the PRE (b).

Figure 3.

The cruise transects and the CTD cast stations: Under northeasterly/easterly winds (a), under easterly/southeasterly winds (b), under southerly/southwesterly winds (c). The red pentagram shows the location of the Shek Pik tide gauge station.

Table 1 lists the tidal phases and wind conditions during different observed transects. Transects B and D were in flood phases during the spring tide under northeasterly winds, whereas transects C and F were in ebb phases. Transect H was during the flood under the easterly wind. In the neap tide period, transect I was during ebb under the southeasterly wind, whereas Transect J was during on the flood tide. Under the southerly wind, transect M was from flood to ebb, whereas transects N and O were during ebb. This data set has high sampling spatial resolution coverage over the PRE in a spring-neap cycle, which is particularly useful to identify the salinity and circulation characteristics and to validate model results.

Date and time (local time)	Transect	Tidal phase	Wind
04 May; 09:43–13:13	Transect B	Flood/spring	Northeasterly
04 May; 13:53–16:16	Transect C	Ebb/spring	Northeasterly
05 May; 09:01–11:51	Transect D	Flood/spring	Northeasterly
05 May; 14:22–16:45	Transect F	Ebb/spring	Northeasterly
07 May; 07:22–15:09	Transect H	Flood/neap	Easterly
08 May; 08:58–12:38	Transect I	Ebb/neap	Southeasterly
08 May; 13:00–15:23	Transect J	Flood/neap	Southeasterly
09 May; 13:50–15:16	Transect M	Flood to ebb/neap	Southerly
10 May; 9:00–11:33	Transect N	Ebb/neap	Southerly
10 May; 12:00–15:01	Transect O	Ebb/neap	Southerly

Table 1.

Observed tidal phases and wind conditions for different transects.

3. Observation results

3.1 Surface salinity

Figure 4 displays the observed surface salinity along the cruise tracks under the wind mainly in the northeasterly during the spring tide from 4 to 7 May (a) and during the neap tide under the wind mainly in the southerly from May 8, 2014 to May 11, 2014 (b). Under the northeasterly wind during the spring tide, freshwater flowed into the coastal sea mainly on the west side of the estuary due to the Coriolis force and the wind forcing, and high salinity water appeared in mid-estuary of the Lantau Channel, implying a high salinity intrusion along the channel. There was a strong salinity gradient in the cross-estuary direction. When the southerly winds dominated during the neap tide, the surface salinity in the estuary is lower than that under the northeasterly wind during the spring tide in the mid-estuary and on the west side, and the surface low salinity water may further spread to the west side. The asymmetry in the surface spatial distribution of the salinity suggests the influence of the wind and tidal forcing on the estuary stratification. In the spring tide and the easterly wind, the low salinity water is confined on the west side and the stratification is weak, while in the neap tide and southerly wind, the stratification is enhanced. Wind forcing is an obvious factor influencing the spatial surface salinity distribution and the stratification [8]. Another important factor is the tidal mixing that is higher in the spring tide and weaker in the neap tide. The stronger turbulent mixing in the spring tide might further decrease the stratification, resulting in higher horizontal salinity gradient in the cross-estuary direction. In the neap tide, the weaker mixing may facilitate the formation of the stratified water in the PRE.

Figure 4.

Observed surface salinity along the cruise tracks under the wind mainly in the northeasterly during spring tide from 4 to 7 May (a) and during neap tide under the wind mainly in the southerly from May 8, 2014 to May 11, 2014 (b).

3.2 Cross-channel salinity distribution and velocity

Figure 5 shows the sectional salinity and current velocity along transects B and C. Transects B and C were under the northeasterly wind during the spring tide on the flood and ebb, respectively. In order to remove high-frequency noises, a 2D rotationally symmetric Gaussian low-pass filter with a (3 × 3) size is used to process the ADCP data. The positive directions of the coordinate system of the current field are in the northward (y) and the eastward (x). Both transects B and C are almost the cross-estuary survey tracks, crossing the Lantau Channel in the mid-estuary. On the flood tide during transect B, the higher salinity water occupied the Lantau Channel, and lower salinity water was located on the two sides ( Figure 5a ). The freshwater was confined on the west side of the estuary; therefore, the salinity on the east side was higher than that on the west side. The isohalines were almost vertical on the west side with relatively strong stratification just in the Lantau Channel. The cross-estuary current was relatively weak with a westward flow at the surface, whereas a velocity divergence appeared in the Lantau Channel, with the eastward (positive) velocity on the east side and westward (negative) velocity on the west side of the Lantau Channel ( Figure 5b ). On the flood tide of transect B, the water flowed into the estuary (northward/landward) with the strongest current on the west side of the Lantau Channel; when the flood tide turned to the ebb tide, the northward current was weaker and reversed on the east side ( Figure 5c ).

Figure 5.

The sectional salinity (practical salinity unit, psu) (a, d), cross-channel circulation (m s⁻¹) (b, e), and along-channel circulation (m s⁻¹) (c, f) along transects B (upper) and C (lower), respectively.

When the tide was turning from flood to ebb during transect C, the salinity reached a maximum (more than 24 psu) in the Lantau Channel ( Figure 5d ) (vs. 21 psu for transect B). Furthermore, the turbulent mixing was enhanced in the Lantau Channel, resulting in weak stratification during transect C. On the ebb in transect C, the freshwater flowed out of the estuary with seaward current in the surface, especially on the west side; at the bottom layer in the Lantau Channel, a landward current appeared for the flood to ebb transition ( Figure 5f ). The cross-channel current was highly complex with alternating convergence and divergence ( Figure 5e ).

During both of the transect B and C periods, the high salinity water flowed into the estuary in the deep Lantau Channel, which could facilitate the appearance of the high-density water in the mid-estuary of the Lantau Channel even during the ebb tide. The survey data indicate that the salinity intrusion could exist on the early ebb for transect C at the bottom of the Lantau Channel during the spring tide and northeasterly wind.

Figure 6 displays the sectional salinity and current velocity along transects D and F. Both transects D and F crossed the Lantau Channel in the lower estuary under the northeasterly wind during the spring tide. In the flood period for transect D, the water flowed into the estuary (northward/landward) with the strongest current appearing in Lantau Channel, reflecting the salinity intrusion along the deep channel, except at the surface on the west side, where the water flowed out of estuary (southward/seaward) influenced by the river discharge ( Figure 6c ). Due to flood tide during transect D, the inflow might bring high salinity water into the estuary especially in the Lantau Channel. This flow pattern might cause the strong density gradient in the cross-estuary direction ( Figure 6a ). The low salinity water resided on the estuary west side, and the high salinity water was located on the east side. In the cross-estuary direction, the westward current appeared throughout the whole depth of the water column, implying that the cross-estuary surface current was induced by the northeasterly wind, while in the deep Lantau Channel, a velocity divergence appeared with eastward (positive) velocity on the east side and westward (negative) velocity on the west side of the deep channel ( Figure 6b ) [10].

Figure 6.

The sectional salinity (psu) (a, d), cross-channel circulation (m s⁻¹) (b, e), and along-channel circulation (m s⁻¹) (c, f) along transects D (upper) and F (lower), respectively.

On the ebb tide during transect F, the maximum seaward (southward) flow in the surface layer reached as large as 0.8 m s⁻¹ with a weak current on the east side of the estuary. However, the bottom water flowed landward (northward) into the estuary, especially in the channel ( Figure 6f ). This was typical of a two-layer structure of density-driven estuarine circulation. Compared with transect D on the flood, the salinity gradient was weak and showed the two-layer structure, although the general pattern of the west low and east high salinity still appeared ( Figure 6d ). The cross-estuary velocity exhibited weak westward flow at the surface and strong eastward flow beneath the surface, related to the wind-driven current in the estuary ( Figure 6e ).

The observation suggests the salinity intrusion existed both for the flood (transect D) and ebb tide (F) same as that for transects B and C. However, different from the observation from transects B and C, the low salinity water further expanded to the east side of the estuary. This may be due to the fact that the seaward flow appeared on the west side both for the flood (transect D) and ebb (transect F) periods.

Transect H surveyed near the estuary mouth and south of Hong Kong Island under the easterly wind. Figure 7 shows the sectional salinity and current velocity along transect H. Near the mouth, low salinity water resided in the whole surface layer (west of 113.87°E) and high salinity water in the bottom layer, while east of 113.87°E, the vertical sectional salinity was well-mixed with the coastal water, indicating that the plume water was appearing most of the surface layer in the estuary mouth ( Figure 7a ). In the along-estuary direction, water flowed into the estuary (northward/landward) during the flood tide ( Figure 7c ). At the surface in the estuary mouth, the water flowed out of estuary (southward/seaward), as it was influenced by the river discharge. In the cross-estuary direction, there appeared a strong westward current on the surface, indicating a wind-driven current induced by the easterly wind. However, a strong eastward current was found in the surface of the estuary mouth due to the plume water spreading, resulting in a velocity convergence at the plume water boundary, indicating a supercritical plume front at the plume boundary ( Figure 7b ).

Figure 7.

The sectional salinity (psu) (a), cross-channel circulation (m s⁻¹) (b), and along-channel circulation (m s⁻¹) (c) along transect H.

During the neap tide period, the wind direction changed from the northeasterly to the southeasterly, and the easterly wind component was weaker than that in the spring tide ( Figure 3b ). The sectional salinity and current velocity along transects I and J are displayed in Figure 8 . During the ebb tide of transect I, the upper layer water flowed out of the estuary (southward/seaward), with a strong current in Lantau Channel, while there was a bottom current flowing into the estuary ( Figure 8c ). In the cross-estuary direction, on the surface of the deep channel, there appeared a strong velocity convergence, an eastward current on the west side and a westward current on the east side, with a weak velocity divergence on the left of the channel at the bottom ( Figure 8b ) [10]. The salinity had a sharp gradient at the velocity convergence location of 113.80°E. The salinity front was more apparent as compared with that during flood tide ( Figure 8a ).

Figure 8.

The sectional salinity (psu) (a, d), cross-channel circulation (m s⁻¹) (b, e), and along-channel circulation (m s⁻¹) (c, f) along transects I (upper) and J (lower), respectively.

In the early flood of transect J, the current was landward/northward (into the estuary) in Lantau Channel, while outside the channel on both sides, the water flowed seaward/southward (out of the estuary), especially in the surface layer ( Figure 8f ). In the cross-estuary direction, a velocity convergence appeared on the west side of the Lantau Channel in the bottom layer ( Figure 8e ). The density structure had still kept an ebb-like structure for a while until flood tide returned in this early flood.

4. Numerical model and settings

4.1 The finite volume community ocean model (FVCOM)

The FVCOM is a 3D unstructured-grid, free-surface, primitive equation, finite volume coastal ocean and estuarine circulation model with triangular grids in the horizontal and terrain-following coordinates in the vertical and wet/dry treatment, developed originally by the University of Massachusetts Dartmouth (UMASSD) and improved by the efforts of Woods Hole Oceanographic Institution (WHOI) [11, 12, 13]. The finite volume approach used in FVCOM, combined with a flexible unstructured grid, provides a better representation of mass, salt, momentum, and heat conservation in an estuarine area with irregular coastline and geography than that used in a rectangular grid system. Mode-split or semi-implicit schemes can choose to solve the governing equations in Cartesian or spherical coordinates. The governing equations, including momentum, continuity, temperature, salinity, and density equations, are closed by using the Mellor and Yamada level 2.5 (MY-2.5) turbulence closure for vertical eddy viscosity [14], or the Smagorinsky eddy parameterization for the computation of the horizontal diffusion [15], as well as an alternative selection of the General Ocean Turbulence Model (GOTM) [16]. For more computational speed, a “mode-splitting” time-stepping scheme is employed to solve the integration process, dividing into internal and external modes with two different time steps. FVCOM has been widely used in estuarine circulation and river plume dynamic studies worldwide, such as the Pearl River Estuary [17], Changjiang River [18, 19], Okatee River [20], Puget Sound estuarine system [21], and Tampa Bay [22].

4.2 FVCOM settings

In the FVCOM, the horizontal grid uses unstructured triangular cells, and the realistic topography is represented using terrain-following coordinates. The greatest advantage of this model is its geometric flexibility since it uses the unstructured and triangular grid meshes, which can provide a real fitting of the irregular coastal boundary. Figure 9a shows the unstructured triangular model grids of the FVCOM for the PRE and coastal water regions. The model domain (111.5–116.5°E, 20–23°N) covers the entire Pearl River Estuary with an open boundary in the northern South China Sea. The total node number in the model grid is 41,027 with 78,539 cells. The horizontal grids have spatial resolutions that vary from 0.1 km to 10 km over the entire domain, with 0.1–0.3 km inside the Pearl River Estuary, 0.3–0.5 km in the estuary mouth, 1.0–2.0 km in the Guangdong coastal water, and 10 km close to the open boundary. The vertical coordinate has 20 levels in uniform hybrid terrain-following grids. The numerical simulations are run for the period from March 1, 2014 to May 30, 2014, during which the observations obtained from the cruise survey conducted between May 3, 2014 and May 11, 2014 were available for model comparisons.

Figure 9.

Model domains for the FVCOM (a) and the WRF model (b).

Tidal forcing is applied to the open boundary including the eight major tidal constituents of M2, N2, S2, K2, K1, O1, P1, and Q1, which are obtained through interpolation from the 1/6° inverse tidal model results of [23], together with the salinity, temperature, and velocity on the open boundary from Hybrid Coordinate Ocean Model (HYCOM) (https://hycom.org/dataserver/glb-analysis) outputs. The meteorological forcing including the wind stress, net heat flux, and evaporation-precipitation (E-P) balance is generated by the Weather Research and Forecasting (WRF) Model. The monthly means of freshwater discharges from the eight river inlets constitute the lateral boundary conditions on the land side [2, 7]. The model is initialized with the March climatology of salinity and temperature fields derived from the World Ocean Atlas 2009 (WOA2009) (https://www.nodc.noaa.gov/OC5/WOA09/pr_woa09.html) and spun up from zero velocity and an undisturbed sea surface elevation. The topographic data for the FVCOM model are derived by interpolation of an electronic navigation chart data for the estuary and coastal water area, and for the offshore region, they are obtained from the global bathymetry data of the General Bathymetric Chart of the Oceans (GEBCO) (http://www.gebco.net/) with 30 seconds (1/120°) horizontal resolution.

The temporal integration of the model uses a split-mode time-stepping method with a 2-second external time step and a split number of 5. The wet/dry treatment is used since some upstream coastal areas can be inundated in high tide. The MY-2.5 and Smagorinsky turbulent closure schemes are used for vertical and horizontal mixing, respectively [14, 15].

In addition to the ocean model, an atmospheric model, the WRF Model, is implemented to provide high spatial and temporal resolution atmospheric forcing for the FVCOM. The WRF Model is supported and developed by the National Center for Atmospheric Research (NCAR), the National Centers for Environmental Prediction (NCEP), the Air Force Weather Agency (AFWA), the Naval Research Laboratory, the University of Oklahoma, and the Federal Aviation Administration (FAA) in the USA. This model is a next-generation mesoscale numerical weather forecasting system that was designed to serve both weather forecasting and atmospheric research needs [24]. A dynamical downscaling technique is employed in the WRF for the regional atmospheric modeling to obtain a higher-resolution output. To realize this, a three-domain nested configuration is designed, as shown in Figure 9b . The outer domain, d01 for the atmospheric model, covers the western Pacific Ocean, the entire South China Sea, and the eastern Indian Ocean, with a horizontal resolution of 72 km; the middle domain, d02, covers the entire South China Sea and southern China, with a horizontal resolution of 24 km; the inner domain, d03, covers the northern South China Sea with a horizontal resolution of 8 km. All domains have 27 layers in the vertical. The NCEP FNL Operational Global Analysis data (http://rda.ucar.edu/datasets/ds083.2/) with a horizontal resolution of 1° × 1° is used to provide initial conditions and lateral boundary conditions for the outer domain. The physics options of the WRF model include the Ferrier microphysics scheme [25], RRTM longwave radiation scheme [26], Dudhia shortwave scheme [27], YSU PBL scheme [28], and Kain-Fritsch cumulus scheme [29]. The model simulation period is from April 20, 2014 to May 30, 2014 covering the cruise survey period. The sea surface winds from the inner domain of WRF model are used to drive the ocean model.

5. Model validation

For the modeling, the coordinate is shown in Figure 9a with the eastward as the positive x-axis and u velocity, the northward as the positive y-axis and v velocity, and the upward as the positive z-axis and w velocity. Since the PRE is generally along the meridional direction, we choose the u to represent the cross-estuary velocity and the v the along-estuary velocity.

The comparison of the cruise observations and FVCOM modeling salinity profiles at a number of stations (B3, D4, D6, F2, G2, G3, H3, H8, I5, J3, M3, and N7 illustrated in Figure 2 ) is shown in Figure 10 . It reveals good agreements between the observations and modeled salinity profiles except at B3, D6, H8, J3, and N7, where the model underestimates the salinity profiles. B3, D6, and N7 were close to coastal land (or island), which may affect the simulation results. Another possible reason is that the river discharge used for the modeling is the climatological data that may not represent the actual discharges, which can cause biases in the salinity profiles between the modeling results and the observations.

Figure 10.

Comparison of the modeling salinity profiles (blue lines) with the cruise observation (red lines).

The comparison of sea levels between the model results and observations at eight tide gauge stations from 1 to 30 May is shown in Figure 11 . The locations of the tide gauge stations are shown in Figure 12 . The hourly tide gauge sea level data during the period between May 1, 2014 and May 3, 2014 are provided by the Hong Kong Observatory (HKO) and Marine Department for the model validation. The modeled sea levels agree well with the observed records both at spring and neap tides for this period.

Figure 11.

Comparison of modeled (blue lines) and observed (red lines) sea surface elevations at eight tidal gauge stations.

Figure 12.

The locations of the tide gauge stations.

To further quantify the accuracy of model results, four statistical parameters are calculated, namely, root-mean-square (RMS) error, relative average error (RE), correlation coefficient, and the model skill parameter used by [30].

The RMS error is defined as:

$$\mathrm{RMS}={\left\{\frac{1}{\mathrm{N}}\sum _{\mathrm{i}=1}^{\mathrm{N}}{\left({\eta }_{\mathrm{mo}}‐{\eta }_{\mathrm{ob}}\right)}^2\right\}}^{1/2}$$ E1

The relative average error (E) is defined as:

$$\mathrm{E}=\frac{{\sum }_{\mathrm{i}=1}^{\mathrm{N}}{\left({\eta }_{\mathrm{mo}}‐{\eta }_{\mathrm{ob}}\right)}^2}{{\sum }_{\mathrm{i}=1}^{\mathrm{N}}\left\{{\left({\eta }_{\mathrm{mo}}‐{\overline{\eta }}_{\mathrm{mo}}\right)}^2+{\left({\eta }_{\mathrm{ob}}‐{\overline{\eta }}_{\mathrm{ob}}\right)}^2\right\}}\times 100\%$$ E2

The correlation coefficient is defined as:

$$\mathrm{R}=\frac{{\sum }_{\mathrm{i}=1}^{\mathrm{N}}\left\{\left({\eta }_{\mathrm{mo}}‐{\overline{\eta }}_{\mathrm{mo}}\right)\left({\eta }_{\mathrm{ob}}‐{\overline{\eta }}_{\mathrm{ob}}\right)\right\}}{{\left[{\sum }_{\mathrm{i}=1}^{\mathrm{N}}{\left({\eta }_{\mathrm{mo}}‐{\overline{\eta }}_{\mathrm{mo}}\right)}^2{\sum }_{\mathrm{i}=1}^{\mathrm{N}}{\left({\eta }_{\mathrm{ob}}‐{\overline{\eta }}_{\mathrm{ob}}\right)}^2\right]}^{1/2}}$$ E3

The model skill parameter is:

$$\text{Skill}=1-\frac{{\sum }_{\mathrm{i}=1}^{\mathrm{N}}{\left({\eta }_{\mathrm{mo}}‐{\eta }_{\mathrm{ob}}\right)}^2}{{\sum }_{\mathrm{i}=1}^{\mathrm{N}}{\left(\left|{\eta }_{\mathrm{mo}}‐{\overline{\eta }}_{\mathrm{ob}}\right|+\left|{\eta }_{\mathrm{ob}}‐{\overline{\eta }}_{\mathrm{ob}}\right|\right)}^2},$$ E4

where η _mo and η _ob represent the model sea level data and the observations, respectively, and the overbar denotes the temporal average. N is the number of records.

Table 2 provides the validation results at the eight tidal gauge stations. The RMS errors are less than 0.15 m, and the relative average errors are less than 3% except at the Ko Lau Wan Station. The correlation coefficients are higher than 0.97, and the skill parameters reach 0.98. Based on the skill assessment results, the model works well in simulating the sea level variations.

	Kwai Chung	Cheung Chau	Shek Pik	Ma Wan	Ko Lau Wan	Tai Miu Wan	Waglan Island	Quarry Bay
RMS (m)	0.11	0.09	0.10	0.09	0.15	0.10	0.11	0.10
RE (%)	2.8	2.1	2.1	2.0	5.3	3.5	2.9	2.9
CC	0.97	0.98	0.99	0.98	0.95	0.97	0.97	0.97
Skill	0.99	0.99	0.99	0.99	0.97	0.98	0.99	0.99

Table 2.

The RMS error, RE, correlation coefficient (CC), and the skill assessment parameter (skill).

Figure 13 presents the comparison of cross-estuary and along-estuary velocities between the cruise observations and FVCOM modeling for transects D and I. The observation duration of transect D was in the flood tide with the northeasterly wind. In the cross-estuary direction, for transect D, the modeled velocity agrees well with the observation, which reveals that the lateral current flows westward (negative) on the left side of the Lantau Channel and eastward (positive) on the right, forming a divergence on the deep channel, although the modeled current speed has some biases from the observed speed, especially for the eastward flow in the Lantau Channel where the modeling underestimates the current ( Figure 13a vs. Figure 13b ). In the along-estuary direction for transect D, the modeled current flows northward with maximum velocity in the surface at the Lantau Channel location, in accordance with the observed current ( Figure 13c vs. Figure 13d ). For transect I with the southeasterly wind during the ebb tide, there was a velocity convergence occurred in the Lantau Channel in the cross-estuary direction as revealed in both observations and modeling ( Figure 13e vs. Figure 13f ). In the along-estuary direction, the modeled velocity agrees well with the observations, and both of the observations and modeling show that the exchange flows exist along the estuary with the strong seaward current in the surface and weak landward current in the bottom ( Figure 13g vs. Figure 13h ). Since the smoothed topography is applied to the modeling, the simulated velocity is spatially less varied than that from the observation; however, the general patterns of the modeling velocity and the observed velocity are similar and consistent, indicating that the FVCOM modeling can well simulate the estuarine water dynamics of the PRE.

Figure 13.

The cruise-observed cross-estuary (a for Transect D, e for Transect I) and along-estuary (c for Transect D, g for Transect I) velocities vs. the model cross-estuary (b for Transect D, f for Transect I) and along-estuary (d for Transect D, h for Transect I) velocities.

6. Variation of estuarine circulation in the PRE

In order to identify the variation of the longitudinal circulation in the spring and neap tide cycle in the PRE, the wind forcing is withdrawn in a modeling run of the FVCOM in the PRE, and modeling results without wind forcing are generated. To remove tidal oscillations, a 25-hour low-pass filter is applied to the modeling outputs of the velocity, salinity, and eddy viscosity. An along-estuary and cross-estuary sections are selected to display modeling results, as shown in Figure 14 .

Figure 14.

The Pearl River estuary with bathymetry contours (m). The along-estuary and cross-estuary sections are shown in blue lines. The intersection is marked as A.

6.1 Along-estuary circulations

Figure 15 illustrates the tidally averaged along-estuary velocity, salinity, and vertical viscosity for the spring and neap tides along the Lantau Channel section. A prominent two-layer exchange flow is shown, with seaward current in the surface and landward current in the bottom ( Figure 15a , b ). The maximum surface seaward current speed is 0.37 m s⁻¹ during the neap tide about 42% stronger than that (0.26 m s⁻¹) during the spring tide in the lower estuary. The bottom landward current is weaker in the spring tide than that in the neap tide. This suggests that the exchange circulation is stronger in the neap tide than that in the spring tide. However, although the bottom landward current is weaker in the spring tide, it may extend farther landward toward the estuary head ( Figure 15a vs. Figure 15b ).

Figure 15.

Tidally averaged along-estuary circulation (m s⁻¹) (a, b), salinity (pus) (c, d), logarithm of vertical viscosity [log(m² s⁻¹)] (e, f) during spring tide (left) and neap tide (right) without wind forcing.

Vertical salinity difference of 2–4 psu appears in the estuary during the spring, whereas it is 4–10 psu during the neap tide, revealing that the stratification during the neap tide is stronger than that during the spring ( Figure 15c vs. Figure 15d ). The maximum stratification is located in the middle of the estuary, and the well-mixed freshwater appears in the upper estuary from 22.5°N during both the spring and the neap tide. It is obvious that the surface salinity is higher in the spring tide than that during the neap tide; however, in the bottom, the salinity intrusion during the neap tide is stronger than that during the spring tide. Without the wind forcing, the estuarine exchange flow pattern is more controlled by the mixing inside the estuary.

As shown in Figure 15e and f , the bottom vertical viscosity is higher than that on the surface both during the spring and neap tides, suggesting that the mixing is generated by the bottom friction that applies to the tidal current inside the PRE [6]. The vertical viscosity is higher during the spring tide than that during the neap tide, indicating that the turbulent mixing is more energetic during the spring tide. The averaged viscosity over the along-channel section is 0.0041 m² s⁻¹ during the spring and is 0.0026 m² s⁻¹ during the neap tide. The PRE shows a classical pattern of the spring-neap cycle of the along-estuary circulation which is mainly controlled by the mixing inside the estuary. The higher mixing dissipates more kinetic energy of the residual current, resulting in a reduced exchange flow in the PRE. Although the horizontal density gradient is stronger during the spring tide, the vertical diffusion may overwhelm the forcing aroused by the density gradient, leading to the weaker longitudinal circulation on the spring tide than on the neap tide.

6.2 Cross-estuary circulations

Figure 16 shows the tidally averaged cross-estuary velocity, salinity, and vertical viscosity for the spring and neap tides across the Lantau Channel. The freshwater outlets are mainly located on the west side of the PRE, which helps to build up a density gradient in the cross-estuary direction, and therefore, a clockwise gravitational circulation (eastward current in the surface and westward current in the bottom) can be formed west of the Lantau Channel as shown in Figure 16a,b (looking into the estuary). The clockwise circulation is stronger in the neap tide than that in the spring tide. Due to the weaker westward current on the east side of the estuary, the eastward current in the surface may extend farther eastward in the spring tide ( Figure 16a vs. Figure 16b ). On the other hand, during the neap tide, there is a convergence in the surface layer located on the east side of Lantau Island, and the stratification is weak on the west side of estuary due to that the well-mixed freshwater reduces the vertical salinity difference ( Figure 16c , d ). However, in the Lantau Channel, the vertical salinity difference is greater in the neap tide, revealing that the stratification during the neap tide is stronger than that during the spring ( Figure 16c vs. Figure 16d ). The vertical viscosity is higher during the spring tide than that during the neap tide, especially on the west side of estuary with shallow water, indicating that the turbulent mixing is more energetic during the spring tide and the circulation is controlled by the tidal mixing over the shoal ( Figure 16e vs. Figure 16f ).

Figure 16.

Tidally averaged cross-estuary circulation (m s⁻¹) (a, b), salinity (pus) (c, d), logarithm of vertical viscosity [log(m² s⁻¹)] (e, f) during spring (left) and neap (right) tide without wind forcing.

7. Discussions and conclusions

The investigation on the Pearl River estuarine circulation is pursued toward better understanding the Pearl River estuarine circulation dynamics. A hydrographic cruise survey for collecting in situ data of marine environmental parameters in the PRE was conducted in the period between May 3, 2014 and May 11, 2014.

The observation shows that the salinity intrusion along the Lantau Channel exists both for the flood and ebb. The intruded salinity should be responsible for the appearance of the high surface salinity water in the channel located inside the estuary, which facilitates the formation of the “sandwich” structure of the horizontal salinity distribution (high salinity in the channel location and low salinity on the east and west sides). However, instead of the “sandwich” type of the horizontal salinity distribution, the lower estuary shows a two-layer structure of the salinity in the horizontal cross-estuary direction, which may be associated with the stronger salinity intrusion in the lower estuary. The observation reveals that there is an exchange flow structure in the longitudinal direction in the PRE with the surface water flowing out of the estuary and bottom water flowing into the estuary, especially along the Lantau Channel. Another important observation is that there is longitudinal convergence or divergence of the cross-estuary velocity close to the channel location for certain tidal conditions. The observation also indicates that the circulation and density structures are highly varied due to the difference in the wind and tidal forcing for different periods.

The numerical simulations of the FVCOM in the PRE and the adjacent water area are conducted with a horizontal spatial resolution as high as 100 m inside the estuary and 20 vertical levels. An atmospheric model, the WRF Model is also implemented to provide with high spatial and temporal resolution atmospheric forcing for the FVCOM. The comparison of model outputs with in situ observations of salinity and velocity profiles and sea level suggests that the modeling in the PRE can well simulate the estuarine circulations. In a model experiment with tide excluding wind forcing, the PRE longitudinal exchange circulation exhibits a strong spring-neap cycle with the maximum circulation in the neap tide and the minimum in the spring tide, revealing that the mixing induced by the tidal current is a dominant factor influencing the circulation. The higher mixing dissipates more kinetic energy of the residual current, resulting in a reduced exchange flow in the PRE.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (grant number 41376035), the General Research Fund of Hong Kong Research Grants Council (RGC) (grant numbers CUHK 14303818, 402912, and 403113), and the talent startup fund of Jiangxi Normal University.