A Numerical Simulation of the Seasonal Cycle of Temperature, Salinity and Velocity Fields in Tokyo Bay

A Numerical Simulation of the Seasonal Cycle of Temperature, Salinity and Velocity Fields in Tokyo Bay

PII: Marine Pollution Bulletin Vol. 43, Nos. 7±12, pp. 145±153, 2001 Ó 2001 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0025-3...

322KB Sizes 2 Downloads 94 Views

PII:

Marine Pollution Bulletin Vol. 43, Nos. 7±12, pp. 145±153, 2001 Ó 2001 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0025-326X/01 $ - see front matter S0025-326X(01)00116-3

Reports

A Numerical Simulation of the Seasonal Cycle of Temperature, Salinity and Velocity Fields in Tokyo Bay FUMIO HORIGUCHI *, JOJI YAMAMOTOà and KISABURO NAKATA§  Research Center for Chemicals Risk Management, National Institute of Advanced Industrial Science and Technology, 16-1 Onogawa, Tsukuba, 305 8569 Ibaraki, Japan àFuyo Ocean Development and Engineering Co. Ltd., 3-15-7 Kuramae, Taitou, Tokyo, Japan §University of Tokai, 3-20-1 Orido, Shimizu, Shizuoka, Japan To understand when oxygen-depleted waters occur, how they develop and when they dissipate in inner Tokyo Bay, realistic simulations were attempted with ®ne spatial and temporal resolution by applying realistic time dependent external forcing. A 3D hydrodynamic model was driven by time-dependent external forcing factors/parameters such as solar radiation, air temperature, relative humidity, wind velocity, and ¯uvial discharge, under the open boundary conditions of 1995. A simulated time series of salinity and temperature agreed fairly well with observed data, except in summer. The model failed to reproduce the development of the surface mixed layer in summer. Several sensitivity analyses on the external forcing parameters such as wind velocity and vertical di€usivity were conducted to reproduce the mixed layer. However, changing these parameter values did not improve the model results. Ó 2001 Elsevier Science Ltd. All rights reserved. Keywords: coastal zone; currents; hydrodynamics; simulation; physical ®elds; Tokyo Bay.

Introduction In Japan, chemical oxygen demand (COD) concentration has been adopted as an environmental standard in coastal waters. Governmental controls for COD emissions from industrial and sewage treatment facilities have been implemented at ®ve year intervals since 1978. However, the water quality of Tokyo Bay, which is

*Corresponding author. Tel.: +81-298-61-8378; fax: +81-298-618357. E-mail address: [email protected] (F. Horiguchi).

typical of semi-enclosed, eutrophied bays in Japan, still shows higher COD concentration than those allowed by the environmental standard. The inner area of Tokyo Bay receives a large amount of inorganic and organic pollutants from the many industrial facilities and huge urban area in the accompanying drainage basin. It has been suggested that the oxygen-depleted waters that are observed in Tokyo Bay during summer are associated with strati®cation of the water column. It is well known that low oxygen or anoxic conditions in bottom waters combined with the increase of sulfur compounds in the bottom sediments can cause serious damage to the benthic fauna. Additionally, under certain environmental conditions, upwelling of these oxygen-depleted waters occur along shoreline, which results in waters overlying tidal ¯at areas and shallower coastal regions in having lower oxygen concentrations. This well known phenomena has been termed `Aoshio', which means the formation of blue color surface water, is due to the in¯uence of sulfur molecules in the low oxygen water. The occurrence of `Aoshio' results in considerable mortality of the short necked clam, which is the most important shell®sh ®shery in this area. From the viewpoint of benthic ®sheries and ecosystem, dissolved oxygen (DO) concentration is a more important environmental health measure than COD. Nakata and Kuramoto (1992) analysed the processes of the formation of oxygen-depleted water in Tokyo Bay with a coupled hydrodynamic and ecological model. They concluded that the hydrodynamic conditions, such as the vertical di€usivity and current system in the bay are the important primary factors in causing the formation of oxygen-depleted water. They also indicated that `Aoshio' occurred in the inner area of Tokyo Bay 145

Marine Pollution Bulletin

and can be explained by upwelling of an oxygendepleted water mass from the bottom layer, driven by northerly wind. However, their study was based on external forcing averaged over a season, or several months, which is static through time. Horiguchi and Nakata (1995) investigated the water quality and carbon cycle in the ecosystem of Tokyo Bay with seasonally averaged meteorological data and principal tidal forcing. They were able to reproduce the spatial distribution of water quality parameters such as COD and nutrients in Tokyo Bay by using a hydrodynamic±ecological coupled model. By using a hydrodynamic±ecological coupled model, the basic mechanism of how the oxygen-depleted water is formed, or how `Aoshio' appears along the coast, has been determined. Understanding when oxygen depletion occurs, how it develops, and when it dissipates in inner-Tokyo Bay, requires more realistic simulations with ®ner spatial and temporal resolution by applying realistic time dependent external forcing. Taguchi et al. (1999a,b) attempted a numerical simulation to analyse seasonal variation in physical and ecological processes in Ise Bay by applying daily based, realistic, time-dependent, external-forcing for three years and succeeded in explaining the history of oxygen-depleted water masses. Nakata et al. (2000) applied the same type of model as Taguchi et al. (1999a,b) to Lake Shinji and Nakaumi. However, there is no such time- dependent simulation study for one year in the case of Tokyo Bay. The present study is concerned with the time dependent structure in the temperature, salinity and velocity ®elds in Tokyo Bay driven by meteorological, tidal and ¯uvial discharge data observed in 1995 as a ®rst step of the analysis.

Governing Equation The governing equations are as follows: 1. Equation of motion ov ov ‡ …v  r†v ‡ w ‡ f0 k  v ot oz Z 0 g ˆ grf rq dz ‡ ‰r  …AH r†Šv q0 H   o ov AZ ‡ : oz oz

…1†

2. Equation of continuity rv‡

ow ˆ 0: oz

…2†

3. The condition at the free surface is given by the equation, Z f  of ‡r vdz ˆ 0: …3† ot H 146

4. The equation for heat conservation is oT oT ‡ …v  r†T ‡ w ˆ ‰r  …KH r†ŠT ot oz   o oT KZ ‡ : oz oz

…4†

5. The equation for salinity conservation is oS oS ‡ …v  r†S ‡ w ˆ ‰r  …DH r†ŠS ot oz   o oS DZ ‡ : oz oz

…5†

6. Equation of state for seawater rt ‡ 1; 1000 X rt ˆ ‡…r0 ‡ 0:1324†f1 qˆ

At ‡ Bt …r0

0:1324†g;

t

r0 ˆ 0:069 ‡ 1:4708S 0:00157S 2 ‡ 0:0000398S 3 ; 2 X …T 3:98† T ‡ 28:30 ;  ˆ 503:570 T ‡ 67:26 t At ˆ T …4:7869 Bt ˆ T …18:03

0:098185T ‡ 0:0010843T 2 †  10 3 ; 0:8164T ‡ 0:01667T 2 †  10 6 ;

…6†

where v is the horizontal velocity, w the vertical velocity, f0 the Coriolis parameter, q the water density, f the sea surface height from the mean sea level, H the depth from the mean sea level to the bottom, T temperature, S salinity, g gravitational acceleration, AH and AZ are the coecients of horizontal and vertical di€usion for momentum, KH and KZ are the coecients of horizontal and vertical di€usion for heat, DH and DZ are the coecients of horizontal and vertical di€usion for salinity, t the time, z the vertical coordinate, positive upward, k the unit vector in the vertical, r is the two-dimensional model operator and q0 density which does not depend on the temperature and salinity with the Boussinesq approximation. The boundary conditions are written as follows. At the surface: 1 S qa 2 q s ˆ ca Wx Wx2 ‡ Wy2 ; q x q …7† 1 S qa 2 q sy ˆ ca Wy Wx2 ‡ Wy2 ; q q oT ˆ Q; …8† qcj oz where sS is the wind stress vector, qa is the air density, c2a is the drag coecient, Wx; y is the velocity of the wind vector of x; y direction, Q is the downward heat ¯ux, c is the speci®c heat of seawater, j is karman constant, F is the fresh water ¯ux (de®ned by the precipitation minus evaporation) and SS is the surface salinity.

Volume 43/Numbers 7±12/July±December 2001

At the bottom: q 1 B sx ˆ c2B uk u2k ‡ v2k ; q q 1 B sy ˆ c2B vk u2k ‡ v2k ; q

…9†

where s2x;y is the bottom friction and uk ; vk is the x; y direction current velocity of bottom level(k). There is neither a heat ¯ux nor a salinity ¯ux across the bottom. In coastal regions, vertical eddy di€usivity and eddy viscosity become important factors. In our model, the following 1-equation turbulent closure model has been employed to evaluate these coecients (Nakata et al., 1983a,b), o o o …hk ek † ˆ …Mk ek † …Nk ek † …xe† zˆ Hk 1 ot ox oy   o oek ‡ …xe† zˆ H ‡ hk D x k ox ox     o oek oe hk D y ‡ Ee oy oz zˆ Hk 1 oy   oe ‡ Ee ‡ Sek Dke ; …10† oz zˆ Hk where

,  2 g oq ou Ri ˆ ; q oz oz p e ˆ L e; p L ˆ k 0 z 1 z=H ;   p g L2 oq E ˆ qL e  exp m ; q e oz

where Ri is the Richardson number, e is the vertical momentum transfer coecient, e is the sub-grid scale (SGS) energy, L is the scale of the length, j0 is the Von Karman constant, z is the distance from the bottom and H is depth, E is the momentum vertical exchange coef®cient, m is the Mamayev constant (m ˆ 0:75 is selected in this study), Dx , Dy are the horizontal dispersion coecient, Ee is vertical exchange coecient and Sek and Dke are the generation and dissipation term of SGS energy in level-k.

Model Run The model bay is shown in Fig. 1. The grid size is 1 km. In Table 1, the computational conditions of the hydrodynamic model are shown. External forcing data are as follows: Meteorological data are daily solar radiation, daily cloud amount, daily air temperature, daily relative humidity, and hourly wind data recorded at an inland weather station in Chiba as shown in Fig. 1 (Japan Meteorological Agency, 1995a). These data are commonly used at all the grid points.

Fig. 1 Model bay. 16 rivers (H), 2 temperature and salinity stations (), an inland weather station (j), and a tidal station ( ). TABLE 1 Computational conditions of the hydrodynamic model. Contents

Selected values

Area

Whole area of Tokyo bay (see Fig. 1) Number of horizontal 72 (north±south)46 (east±west) meshes (see Fig. 1) Grid interval 1 km Layer locations Number of vertical level is 10 Level 1 surface to )2 m, Level 6±10 to )13 m Level 2±2 to )4 m, Level 7±13 to )16 m Level 3±4 to )6 m, Level 8±16 to )19 m Level 4±6 to )8 m, Level 9±19 to )22 m Level 5±8 to )10 m, Level 10±22 to bottom Sea level Tide at Mera (see Fig. 1) Wind condition at Chiba weather station (see Fig. 1) East±west comp.: 9±15 m s 1 North±south comp.: 15±20 m s 1 River conditions 16 rivers (see Fig. 1) Temperature: 9±30°C (observed values) Salinity: 0.0 PSU (observed values) Meteorological Solar radiation, cloud amount, temperature, conditions relative humidity (see Fig. 2) Model parameters Coriolis parameter 8:610 5 s 1 Horizontal eddy vis30 m2 s 1 cosity and di€usivity Vertical eddy viscosity 30 m2 s 1 and di€usivity Equation of SGS Calculated period 365 days

Fluvial water data are daily ¯ow rates of the 16 rivers (Ministry of Construction, 1995). The ¯uvial water salinity is taken as 0 and daily temperature data are used. 147

Marine Pollution Bulletin

Sea level data at the open boundary are hourly data that were recorded at Mera tidal station outside the bay. The location of Mera is shown in Fig. 1 (Japan Meteorological Agency, 1995b). The amplitude of M2 tidal constituent is ca 50 cm, and therefore the strong tidal mixing cannot be expected in Tokyo Bay. Temperature and salinity data at the open boundary are semi-monthly data collected near the model open boundary (Kanagawa Prefectural Fisheries Research Institute, 1995).

The surface heat ¯ux is comprised of solar radiation, long wave radiation, sensible and latent heat ¯uxes, each of which is calculated with the meteorological data and the predicted surface water temperature. The heat ¯ux is assumed to be completely absorbed in the uppermost layer. The downward long wave radiation is calculated with a formula by Brunt. The wind stress is calculated with a bulk formula by using the wind velocity data. Fig. 2 shows the daily solar radiation, cloud amount, air temperature and relative humidity during 1995.

Fig. 2 Solar radiation, cloud cover, air temperature and relative humidity in 1995.

Fig. 3 River discharge rate, shown only for 3 of the 16 rivers: Arakawa, Edogawa and Sumidagawa Rivers.

148

Volume 43/Numbers 7±12/July±December 2001

Fig. 3 shows the river discharge rate for 3 of the 16 rivers. In 1995, the rainy season started in the middle of June and terminated at the end of July, which was re¯ected by the ¯ow rates of rivers in Fig. 3. The e€ect of a typhoon was also seen in the ¯ow rate data on September 20 in Fig. 3. The external forcing is updated at each time step with a linear interpolation in time. The time step is 10 s.

The coecient of vertical di€usion is calculated in the model with an SGS equation. The Coriolis parameter is taken as 8:6  10 5 s 1 and the b e€ect is neglected.

Results and Discussion The model was run with the external forcing on 1 January 1995 after a spin-up calculation of the model. The

Fig. 4 Simulated and observed temperature and salinity at Stn. 01. (a) ®rst layer (surface to 2 m), (b) second layer (2±4 m), (c) ®fth layer (8±10 m) and (d) sixth layer (10±13 m).

149

Marine Pollution Bulletin

time integration is forwarded until the end of the year 1995 to get an annual cycle of the temperature, salinity and velocity ®elds. The output is stored at 1 h intervals.

The simulated temperature and salinity values were compared with observed ones at two stations viz: Stn. 01 and Stn. 02 (shown in Fig. 1) in Figs. 4 and 5. At Stn. 01,

Fig. 5 Simulated and observed temperature and salinity at the ®rst layer (surface to 2 m) at Stn. 02.

Fig. 6 Vertical distributions of simulated and observed temperature and salinity at Stn. 01. (a) August 6, (b) August 9, and (c) August 10. The simulated data are tentatively assigned at middepths of the layers.

150

Volume 43/Numbers 7±12/July±December 2001

the observed data were available for almost every day at every 0.5 m depth between June 1995 and April 1996, while at Stn. 02 they were available only at the sea surface or nearby once a month. At Stn. 01, the simulated and observed results agree well except for August and September. In these months, the model predicted lower salinity and temperature than the observations, especially in the surface layer. The model also fails to simulate short period ¯uctuations in these months. At Stn. 02, the model results show good agreement with the observed values. Figs. 6(a)±(c) give vertical distributions of temperature and salinity at Stn. 01 on August 6, 9 and 10. Although the simulated temperature and salinity values are not de®ned at speci®c depths, they are tentatively assigned at mid-depths of the layers. The observed pycnocline deepened by about 5 m over 4 days from August 6 to August 10. This suggests generation of internal waves. However, no distinct pycnocline, or, no development of a surface mixed-layer appears in the simulation. Figs. 7(a)±(d) show correlation coecients between the simulated and observed surface salinity and tem-

perature at Stn. 01 and Stn. 02. Except for the surface salinity at Stn. 02, they are either 0.80 or 0.98. Even for the surface salinity at Stn. 02, if the September data are excluded, the correlation coecient is 0.81 (shown in Fig. 8). Here we discuss why the discrepancies between the model and the observed results appeared. The underestimation of salinity and temperature in August and September can be attributed to the boundary conditions, such as meteorological data and river ¯ow rate data. The river ¯ow rate is usually not observed at the river-mouth, but at 5±20 km further up the river. Therefore, the observed ¯ow rate in the rainy and typhoon season is smaller than the ¯ow rate incoming into the bay. The model ignores the rainfall, though it decreases the surface salinity in the rainy and typhoon seasons. Given this situation, the simulated ¯ow rates were three times greater than the observed ones in the rainy and typhoon seasons. However the model predicted overestimation of salinity even if such ¯ow rates were given. To ®t the model results to the observed ones, more ¯ow rates or direct precipitation on the sea surface may be required.

Fig. 7 Correlation between the simulated and observed data. (a) salinity, (b) temperature at Stn. 01; (c) salinity, (d) temperature at Stn. 02. The solid line shows the regression line.

151

Marine Pollution Bulletin

Fig. 8 Correlation between the simulated and observed surface salinity data at Stn. 02 except in September when the observed salinity drops strikingly. The solid line shows the regression line.

In regard to the meteorological data, the weather station is located inland as mentioned above (Fig. 1). It was assumed that the meteorological data observed there such as cloud cover and wind ®eld, represented the meteorological conditions over the whole Tokyo Bay area. In reality, these conditions show large spatial heterogeneity, but there is no way to evaluate the spatial distribution of cloud cover or the wind ®eld on a daily basis. For the wind ®eld, there is no observation buoy system to evaluate the wind velocity distribution over Tokyo Bay itself. In this study, the wind velocity observed at Chiba was used. These problems may lead to error in the model prediction. To explain why the model failed to predict the development of the surface mixed layer, we attempted two sensitivity analyses. The ®rst examined the change of wind velocity. The model was run with wind velocity at three times what it was in the original run. Fig. 9 shows

Fig. 10 Comparison among three cases with changed Mamayev constant; () indicates salinity at m ˆ 0:75, ( ) indicates temperature at m ˆ 0:75, () indicates salinity at m ˆ 0:0075, (}) indicates temperature at m ˆ 0:0075, (M) indicates salinity at m ˆ 0:15, (±) indicates temperature at m ˆ 0:15.

the comparison between them. The di€erence between the two runs is very small. The vertical di€usivity was calculated by using an SGS equation with a Mamayev constant (in this study). To examine the in¯uence of the Mamayev constant on the vertical pro®le of salinity and temperature, another two runs with di€erent Mamayev constants (0.0075 and 0.15) were attempted. The smaller magnitude of this constant implies the larger di€usivity. The comparisons among the three cases, shown in Fig. 10, suggest that this constant is not e€ective to simulate the development of the surface mixed layer. The temperature and salinity data used as the open boundary conditions were observed at a station away from the open boundary, so their accuracy is poor. In addition, the open boundary condition does distort more or less dynamical processes in the bay. Agreement between the simulated and observed results in Figs. 4 and 5 is satisfactory for the present model even with these boundary conditions.

Conclusion

Fig. 9 Comparison between original run and run with wind velocity at three times; () indicates salinity at original wind, ( ) indicates temperature at original wind, () indicates salinity at three times wind, (}) indicates temperature at three times wind.

152

A numerical model driven by time-dependent external forcing parameters such as solar radiation, air temperature, relative humidity, wind velocity, ¯uvial discharge, under the open boundary conditions of 1995, was attempted to realistically simulate the physical processes in Tokyo Bay. Simulated time series of salinity and temperature pro®le agreed fairly well with observed data, except in the summer months. Therefore, it was possible in this study to realistically simulate the physical processes in Tokyo Bay for the year 1995. The next step to improve the simulation will be to collect data of better quality, at more stations, for the following parameters, which are of primary importance for the model validation: spatial

Volume 43/Numbers 7±12/July±December 2001

distribution of external forcings, temperature, salinity and velocity. The model employed in this study used a 1-equation type turbulent closure model. However, the model failed to reproduce the surface mixed layer. Therefore, we have to examine the k±e type closure model (e.g. Mellor and Yamada, 1982) to reproduce the mixed layer in the future. Following that, the model will be imbedded into a sophisticated ecosystem model. We thank deeply Dr K. Takano, Professor Emeritus of Tsukuba University, for considerable guidance in this study. In addition, Mr Nakane of Science and Technology co. for collected meteorological data used in the model is thanked. Horiguchi, F. and Nakata, K. (1995) Water quality analysis of Tokyo Bay using mathematical model. Journal of Advanced Marine Science and Technology Society 1(1), 71±92 (in Japanese). Japan Meteorological Agency (1995a) Annual Report of Meteorological data (in Japanese). Japan Meteorological Agency (1995b) Tide data of Mera. Tidal Station (in Japanese). Kanagawa Prefectural Fisheries Research Institute (1995) Observation data (Temperature and Salinity) in Tokyo bay (in Japanese).

Mellor, G. L. and Yamada, T. (1982) Development of a turbulence closure model of geophysical ¯uid problems. Review of Geophysics and Space Physics 20, 851±875. Ministry of Construction (1995) Annual Report of Fluvial Discharge. Japan River Association, Japan (in Japanese). Nakata, K., Horiguchi, F., Taguchi, K. and Setoguchi, Y. (1983a) Three-dimensional simulation of tidal current in Oppa Estuary. Bulletin of The National Research Institute for Pollution and Resources 12, 17±36 (in Japanese). Nakata, K., Horiguchi, F., Taguchi, K. and Setoguchi, Y. (1983b) Three-dimensional eco-hydrodynamical model in coastal region. Bulletin of The National Research Institute for Pollution and Resources 13, 119±134 (in Japanese). Nakata, K. and Kuramoto, T. (1992) A model of formation of oxygen depleted waters in Tokyo Bay. Journal of AMTEC 5, 107±132. Nakata, K., Horiguchi, F. and Yamamuro, M. (2000) Model study of Lakes Shinji and Nakaumi ± a coupled coastal lagoon system. Journal of Marine System 26, 145±169. Taguchi K., Nakata, K. and Ichikawa, T. (1999a) A three-dimensional simulation of long-term variability in the ¯ow ®eld and T±S structure in the Ise-Mikawa Bay estuary. Journal of Advanced Marine Science and Technology Society 5(1&2), 37±48. Taguchi K., Nakata, K. and Ichikawa, T. (1999b) A three-dimensional simulation of the lower trophic ecosystem in the Ise-Mikawa Bay estuary using a coupled physical and biochemical model. Journal of Advanced Marine Science and Technology Society 5(1&2), 49±62.

153