Kranji Reservoir is a tropical reservoir located in the north-western region of Singapore. It was formed by damming the estuary of the Kranji River in 1972. The main objective of the project is to develop an integrated three-dimensional hydrodynamics model for simulation and prediction of tropical reservoir flow hydrodynamics.
This report presents the results from the 2-year simulation after assembly and calibration of the Estuary and Lake Computer Model (ELCOM) for Kranji Reservoir. Kranji Reservoir is modeled using 40 x 40m grid resolution and has an estimated spin-up time of 20 days. The residence time is approximately 30.5months, a multiple of the hydraulic detention time which is estimated to be 15months.
Results show that wind is the most important driving force of hydrodynamics. Air temperature and incident solar radiation impact the thermal balance while tributary inflows and outflow have only local impact on the hydrodynamics. Mass balance is crucial in verifying the inputs for ELCOM for valid simulation. Evaporative losses have to be catered for by reducing output flow rate.
In many water supply reservoirs, the most serious water quality problems are those related to the local accumulation of algae. The Public Utilities Board (PUB) is spending an estimated amount of $1.5M annually in reservoir management related to eutrophication avoidance (e.g. via contracted physical removal of aquatic plans, algal scum, chemical dosing, energy cost in aeration etc). Even then recurrent blue-green algae blooms are still found in most months of the year in hypertrophic reservoirs like Kranji, which we are using for our research. [6]
The main objective of the project is to develop an integrated 3D water quality model for simulation and prediction of tropical reservoir flow hydrodynamics. This will be used in assisting management to enhance water quality and in planning water resources; by predicting the impacts of future land use changes or environmental perturbations.
1.2 Objective of Study
In order to track the transport, quantity, quality and fate of water within a reservoir/catchment system, an integrated water quality/quantity model comprising of simulation models will be assembled.
The Reservoir Model will mainly address the problems of eutrophication and dissolved oxygen balance. It can be used to identify and predict the conditions conducive for high suspended solid loads and eutrophication in reservoir waters (i.e. season and catchments characteristics).
The Reservoir Model will be made up of (1) ELCOM (Estuary and Lake Computer Model) - a hydrodynamics model to predict flow velocities and current directions; and (2) CAEDYM (Computational Aquatic Ecosystem Dynamics Model) - a water quality model made up of a column model to predict the flow of nutrients and algal dynamics and a benthic model to predict the exchange of materials between sediments and the water column.
Scope of Work
The 3D ELCOM (Estuary and Lake Computer Model) of Kranji Reservoir is assembled using meteorological data. Calibration of the model is performed to establish the input parameters. This is done by verifying the results with predictions based on prior experience within a structured framework. Several parameters can be analysed from the results obtained from the simulation runs of the ELCOM. Emphasis has been placed on the analysis of temperature variation and tracer concentration for different sections within the reservoir with respect to time.
CHAPTER 2 - REVIEW OF THEORY
& PREVIOUS WORK
2.1 Introduction
Hydrodynamic models simulate the thermodynamics and currents in lakes, estuaries and coastal oceans. They simulate the mixing and transport of heat, sediment, pollutants and biological particles. Hydrodynamic models provide powerful predictive capability that allows the quantitative assessment of both short- and long-term management strategies and environmental changes.
The main focus of this paper is on ELCOM (Estuary and Lake Computer Model) numerical simulation code developed at the Centre of Water Research (CWR) in Western Australia. CWR is recognized as a world-leader in the development and application of one- and three-dimensional hydrodynamic models of aquatic systems ranging from wetlands, lakes and reservoirs to estuaries and coastal ocean.
ELCOM is a three-dimensional hydrodynamics model used for predicting the velocity, temperature and salinity distribution in natural water bodies subjected to external environmental forcing such as wind stress, surface heating or cooling. The three-dimensional hydrodynamic transport and mixing is simulated through solution of unsteady viscous Navier-Stokes equations and scalar transport equations.
Through coupling with CAEDYM (Computational Aquatic Ecosystem Dynamics Model) water quality module, ELCOM can be used to simulate three-dimensional transport and interactions of flow physics, chemistry and biology. The model solves the unsteady Reynolds-averaged equations using a semi-implicit method with quadratic Euler-Lagrange discretization of momentum advection [Casulli and Cheng, 1992] and a conservative ULTIMATE QUICKEST approach for scalar transport [Leonard, 1991].
2.2 Literature Review
Aquatic systems are modelled for the design and assessment of management strategies and formulate intervention strategies for algal blooms, metals, turbidity, colour, suspended sediments, downstream conditions and bio-manipulation. The power of modelling is in prediction, which is based on prior experience, within a structured framework. Upon understanding the flow dynamics and regimes within the water bodies, the dominant processes can be identified and real-time management can be achieved.
There are four main types of models used for water studies for varying objectives and considerations:
OD - mass balance where inflow is balanced with outflow and sedimentation
(e.g. BATHTUB)
1D - horizontally-averaged (e.g. DYRESM-CAEDYM)
The Dynamic Reservoir Simulation Model (DYRESM) is a one-dimensional hydrodynamics model for predicting the vertical distribution of temperature, salinity and density in lakes and reservoirs.
The Computational Aquatic Ecosystem Dynamics Model (CAEDYM) consists of a series of mathematical equations representing the major biogeochemical processes influencing water quality.
DYRESM-CAEDYM couples the one-dimensional hydrodynamics model DYRESM with the aquatic ecological model CAEDYM. This allows for investigations into the relationships between physical, biological and chemical variables in water bodies over seasonal and inter-annual timescales.
2D - laterally-averaged (e.g. CE-QUAL)
3D - full dynamics (e.g. ELCOM-CAEDYM)
The Estuary Lake Computer Model (ELCOM) is a numerical modeling tool that applies hydrodynamic and thermodynamic models to simulate the temporal and spatial behaviour of water bodies with environmental forcing.
ELCOM-CAEDYM couples the three-dimensional hydrodynamics model ELCOM with the aquatic ecological model CAEDYM. This allows for investigations into the spatial and temporal relationships between physical, biological and chemical variables in water bodies, over single events or up to seasonal timescales.
The choice of model used depends on the ultimate objective. Other considerations include time (flood, seasonal etc) and space (horizontal gradients, vertical stratification) scale of interest. It is also dependent on the data available for running and validation and resource availability in terms of computation, time and money.
The water quality and ecological state of an aquatic system depend not only on the transport of sediment, solutes and biological particles but also on the subsequent growth and biogeochemical transformation of those agents. Water quality modelling requires the coupling of a one- or three-dimensional hydrodynamic model with an aquatic ecological model such as CAEDYM. A coupled water quality model such as DYRESM-CAEDYM or ELCOM-CAEDYM can quantitatively simulate the spatial and temporal variation in concentrations of algae, dissolved oxygen, nutrients and other biological agents. The models can be used predict algal blooms and oxygen depletion, and to assess alternative water quality management strategies.
Two advanced computer simulation models will be assembled for the reservoir/catchment system. The hydrology, hydraulics, water quality and sediment transport along the network of channels and water bodies draining the basin and into the reservoir will be simulated. The one dimensional Rainfall-Runoff model will simulate the rate, quantity and quality of surface runoff arising from spatially and temporally varying rainfall input over the watershed basin. The Reservoir Model is to simulate the flow hydrodynamics and water quality in the reservoir. A complete model will be the integration of both models with outputs of the Rainfall-Runoff Model being used as inputs to the Reservoir Model.
The Reservoir Model will mainly address the problems of eutrophication and the dissolved oxygen balance. Eutrophication of reservoirs in Singapore occurs to varying extents, depending on season and catchment characteristics. It is the excessive growth of algae due to primarily significant nutrient loads causing the depletion of dissolved oxygen which is detrimental to aquatic life. The excessive growth of algae can clog filters and reduce filter run hours during treatment. The anaerobic environment can also release certain contaminants and nutrients in the sediments. Control measures like chemical dosing can be used but it is more environment-friendly to manage the inputs.
The Reservoir Model consists of 2 models:
CAEDYM (Computational Aquatic Ecosystem Dynamics Model) is a water quality model. It enables us to have a better understanding of suspended sediment dynamics and nutrient balances of nitrogen, phosphorus and carbon. [Fig 2.2a] Besides that it also allows us to study eutrophication, food-web, aquaculture and pollution. CAEDYM consists of a series of mathematical equations representing the major biogeochemical processes influencing water quality and may be run independently or coupled with hydrodynamic models like DYRESM or ELCOM to simulate the aquatic ecosystem and physical hydrodynamic processes.
The CAEDYM consists of the following modules which can be invoked: Inorganic particles (SSOL1, SSOL2), Oxygen (DO), Organic nutrients (POM, DOM), Inorganic nutrients (NH4, NO3, PO4, SiO2, DIC (and pH)), Bacteria (BAC), Phytoplankton (CHLA/C, internal nutrients, metabolites), Higher biology (Zooplankton, Fish), Benthic biology (Macroalgae, Clams, Macrophytes etc), Pathogens (Cryptosporidium, Coliforms) and Metals (Al, Fe, Mn).
Fig 2.2a Carbon-Oxygen-Nitrogen-Phosphorus Model [9]
The model consists of two components [Fig 2.2b]:
The water column which is subjected to hydrodynamic flows in the vertical and horizontal directions, and the kinetic processes governed by the chemical and biological processes. It consists of interacting systems like the phosphorus, nitrogen and carbon cycle, phytoplankton dynamics, zooplankton dynamics and the dissolved oxygen balance.
The sediments or benthic layer focuses on mineralisation, nutrient dynamics and sediment oxygen demand.
Fig 2.2b Basic Cycling [9]
State variables [Fig 2.2c] will be subjected to physical transport processes (hydrodynamics) so that spatial and temporal concentrations can be predicted through feedback and interaction. The model will form the basis for coupling to the hydrodynamics model.
Fig 2.2c CAEDYM State Variables
(II) ELCOM (Estuary and Lake Computer Model) is used to predict the flow velocities and current directions. It is the three-dimensional hydrodynamics driver for CAEDYM. The model solves the unsteady Reynolds-averaged Navier-Stokes (RANS) equations using a semi-implicit method with quadratic Euler-Lagrange discretization of momentum advection and a conservative ULTIMATE QUICKEST approach for scalar transport. A one-dimensional mixed-layer model is extended to 3D for turbulence closure of vertical Reynolds stress terms.
The ELCOM numerical method takes its basic structure from the TRIM scheme of Casulli and Cheng (1992) with adaptations for accuracy, scalar conservation, numerical diffusion and implementation of a mixed-layer model.
ELCOM is designed as a model of physical processes rather than a simulation of fluid flow. A simulation method will solve the evolution equations to the designed degree of accuracy given a sufficiently fine grid. The accepted test of "sufficient" resolution requires a grid resolution study with at least one numerical simulation performed with twice the resolution as the base grid scale. If the simulation results change dramatically with the grid resolution, then the base grid is considered to be under-resolved. Since maximum practical resolution is a function of computational power, the accepted procedure in numerical simulation is to reduce the Reynolds and/or Froude numbers of the flow under consideration until a well-resolved flow can be simulated with the available computer. The fundamental tenet of flow simulation is that a well-resolved solution to the evolution equations is invariant to further refinement of the grid or change of numerical method.
Three-dimensional (3D) models of geophysical flows with topographic effects are often not suitable for rescaling without distorting the physics under investigation. Thus the attainable resolution is strictly a function of the available computational power. Carefully conducted coarse-grid 3D models are a necessary bridge between highly-resolved three-dimensional simulation of idealized flows (presently obtainable through DNS, LES and RANS for limited Reynolds and Froude numbers), and 1D models of lake dynamics and 2D models of estuaries. The evolution equations derived from first principles (as used in simulations) must be carefully reconsidered in conjunction with the behaviour of the numerical methods applied.
The hydrodynamic model solves the transport equations which are the unsteady Reynolds averaged Navier-Stokes equations and scalar transport equations using Boussinesq and eddy-viscosity approximations, neglecting non-hydrostatic pressure terms. The free surface evolution is governed by an evolution equation developed by a vertical integration of the continuity equation applied to the Reynolds-averaged kinematic boundary condition. Surface thermodynamics are governed by bulk transfer models.
Issues related to stratification are directly addressed with the full incorporation of the vertical distribution of temperature, salinity and density in the solution process. The forcing in the model is due to wind stress, surface heat and mass transfer, light penetration and inflows from the catchment.
2.3 Previous Work
Application of numerical water quality models to evaluate aquatic management is widely used for lakes, reservoirs, rivers, estuaries and coastal zones. Early studies of numerical modelling of lakes were carried out by Simons (1974), simulating wind-induced circulation in Lake Ontario. A 3D finite-difference method and a constant eddy-viscosity model were used to compute the turbulence. The eddy-viscosity was varied until the computational results agreed with the measurements. Good agreement was found for water velocities and water levels.
Wind stresses, surface heating and density currents form the driving energy fluxes of a stratified lake. The basin-scale energy flux from the wind is of particular interest to Imberger et al. (2000) because of its dominant role in setting the thermocline in motion, which, in the absence of inflows and outflows, is the primary energy store for transport and mixing below the wind-mixed layer. Basin-scale internal waves provide the driving forces for vertical and horizontal fluxes in a stratified lake below the wind-mixed layer, thus examined with a hydrostatic, z-coordinate three-dimensional numerical model at coarse grid resolution.
The 3D estuary and lake computer model ELCOM is applied to modelling Lake Kinneret, Israel, and is compared with field data under summer stratification conditions to identify and illustrate the spatial structure of the lowest-mode basin-scale Kelvin and Poincare waves. The model solves the unsteady Reynolds-averaged Navier-Stokes equation using a semi-implicit method similar to the momentum solution in the Tidal, Residual, Intertidal Mudflat (TRIM) code with the addition of quadratic Euler-Lagrange discretization, elimination of vertical diffusion terms in the governing equations and scalar transport using a conservative flux-limited approach. Instead of the conventional eddy-viscosity/diffusivity assumption, a detailed description is provided of turbulence closure for the vertical Reynolds stress terms and vertical turbulent transport using a 3D mixed-layer model parameterized on wind and shear energy fluxes.
Despite certain numerical modelling issues (like the damping of the Kelvin wave, numerical diffusion of the metalimnion(intermediate region of mixing) and sub-grid scale modelling of internal wave effects) that remain unsolved, this approach gives a good representation of the depth of the mixed layer at coarse vertical grid resolutions that allows the internal waves to be energized correctly at the basin scale.
More recent work involving ELCOM was done by Imberger and Hodges (2003). The 3D ELCOM and field data were used to investigate the influence of spatial and temporal variations in wind-forcing on the circulation in Lake Kinneret. Wind blowing over a lake surface form a highly turbulent surface mixing layer. Turbulence rapidly distributes momentum, transferred from wind to water, over the depth of this layer such that (initially) the surface water moves downwind as a slab (Spigel and Imberger, 1980). This work considers the effect of spatial and temporal variability in the wind field on basin-scale motions in lakes.
ELCOM simulations were started from rest with horizontal free surface. The mixing model was developed based on the one-dimensional integral mixing model of Spigel et al. (1986) adapted by way of Hodges et al. (2000) to 3D fixed-grid architecture. Mixing is computed separately in each water column, starting from the free surface and operating downwards through grid layers. In ELCOM, mixing and advection are calculated separately, thus within a computational step, mixing will not limit the growth of shear until after the advection computation.
Internal wave motions are well reproduced by the numerical model when forced with a spatially uniform wind. However, simulated seiche amplitudes are too large and lead observations by 3-10h (for a 24h-period wave) at different locations around the lake. Consideration of spatial variation if wind field improves simulated wave amplitudes, reducing phase error at all stations. This improvement is attributable to a better representation of the horizontally averaged wind stress and can be reproduced with a spatially uniform wind that has the same horizontally averaged wind stress as the spatially varying wind field. However, a spatially varying wind field is essential for simulating mean surface-layer circulation, which is shown to be predominantly directly driven by wind stress curl and topological moment, and hence the horizontal distribution of wind stress is important.
CHAPTER 3 - METHODOLOGY
Kranji Reservoir is located in the north-western region of Singapore, with a catchment area of approximately 54,500,000m2. There are three tributaries inflows (namely Kangkar, Tengah and Pengsiang) leading to the reservoir that has an estimated storage capacity of 15.85 x 106 m3. The reservoir averages about 6.77m in depth with the maximum at 18m.
A time-dependent three-dimensional hydrodynamic model (ELCOM) of Kranji Reservoir was validated and used to understand the physical processes associated with vertical mixing in the reservoir. This work presents a 2-year simulation that uses measured bathymetry, wind velocity, relative humidity, total solar radiation, air temperature, inflow and outflow. The model indicates that wind is the major driving force of horizontal and vertical water movement.
Fig 3a ELCOM-CAEDYM
The 2 year simulation was performed in three stages, starting with a 6-months simulation before restarting ELCOM for subsequent 6-months and 12-months simulation. Meteorological data for 2004 simulation was duplicated from 2003.
Steps for ELCOM simulations [Fig 3a]:
The bathymetry.dat file that provides topography and grid information is prepared after digitizing the map.
Kranji Reservoir is modeled using grid resolution of 40 x 40m in order to consolidate the survey data points from the surface map. The bathymetry file contains 7874 grid cells in the model horizontal grid (127 x-rows and 62 y-columns) and there are 22 depth layers of 1m each. The bathymetry data is digitized through interpolation, rotation and straightening of bathymetry data to fit the main flow along the axis.
The initial bathymetry [Fig 3b] was constructed using grid resolution of 20 x 20m. This proved to be unsuitable as the Courant-Friedrichs-Lewy (CFL) criterion often exceeded the computational limit of 0.75, causing the simulation to terminate.
Fig 3b Initial Bathymetry
Alterations were made to the final bathymetry [Fig 3c] that included increasing the grid size and trimming off the small bends along the inflows. This significantly reduced the simulation run-time to almost half, allowing longer simulations (i.e. 6 months simulation run could be completed in about 25 hours) without compromising on the resolution and results.
Fig 3c Final Bathymetry
Prepare bc.dat file that designates sets of cells for applying boundary conditions
The boundary conditions include labeling and designating the locations of the inflow/outflow area and aeration points on the bathymetry along with the direction of flow (denoted by the face of the cell which inflow flows into).
Run pre-processor (pre_elcom.exe) to produce sparsedata.unf and usedata.unf files after preparing run_pre.dat
Directories for input (bathymetry.day, bc.dat) and output (sparsedata.unf, usedata.unf) files are specified in run_pre.dat for ELCOM to process and deliver.
Configure the simulation run_elcom.dat file
This is the file where configuration controls are set for ELCOM. Time controls include the date to start the simulation, duration of time steps (e.g. 60s time step) and number of time steps (i.e. duration of simulation. e.g. 6 months simulation has an equivalent of 260640 time steps). Input files (located in the 'Infiles' folder) to be used are specified with file extensions for ELCOM to retrieve. Other options with the file are set as default according to CWR.
Specify output in datablock.db file for designating data output operations
The different output groups, group data types and output sets are specified accordingly. Output groups can be denoted as one-dimensional, two-dimensional or three-dimensional [Fig 3d] along with the output intervals. A few possible forms are 1D point profile, 2D curtain (cross sections), sheet (plan view) and general 3D. This work will be concentrating on presenting points (points located at each inflow at outflow and four PUB points in blue) located within the reservoir and curtains (highlighted in 3 colours in Fig 3c, spanning from individual inflow to outflow) that span the length of the reservoir. Group data types refer to the different variables to be simulated (e.g. temperature, tracer and velocity in the x-, y- and z-directions). Lastly, the output sets specify the number of cells within each output groups.
Fig 3d Three-Dimensional Model of Kranji Reservoir
Prepare temporal boundary conditions files(.dat)
Meteorological data containing wind speed/direction, solar radiation, air temperature and relative humidity
Inflow data containing inflow rate and water temperature (of the 3 tributaries)
and initial conditions files (define vertical variation in initial conditions for any transportable scalar using water temperature and salinity).
Kranji reservoir ELCOM was modeled and calibrated using 2003 meteorological data from Changi Meteorological Station as field data for Kranji is unavailable. Inflow data were compiled to obtain a monthly average value (i.e. rainfall x rainfall coefficient x area). Thus this is assumed to be an ideal simulation and there is no observed data to compare with. Results from the simulation can only be verified using predictions based on prior experience within a structured framework.
Run the executable (elcd_sp.exe)
Post-process data (dbconv.exe) and obtain NetCDF files
Visualization using ELCOM-CAEDYM graphical user interface (GUI)
The post-processing tool (dbconv) converts raw output of ELCOM into the standard NetCDF format. The ELCOM-CAEDYM graphical user interface (GUI) plots the standard NetCDF files and the output can be saved as image or movie file.
Results Analysis using Matlab
Matlab can be used to access the NetCDF files in order to generate graphical plots for analysis.
CHAPTER 4 - PRELIMINARY RUNS,
RESULTS & ANALYSIS
4.1 Preliminary Runs
The major external forcings in the simulation period included in this model are wind speed and direction, air temperature, incident solar radiation and tributary inflow/outflow [Fig 4a]. Like many waterbodies in the world, wind is the most important driving force of hydrodynamics. Air temperature and incident solar radiation impact the thermal balance and tributary inflows and outflow have only local impact on the hydrodynamics.
Numerous simulations were carried out to calibrate and verify ELCOM for the 2-year simulation with tracer input of 30units concentration in Kangkar tributary. The initial problem was exceeding the Courant-Friedrichs-Lewy (CFL) criterion. In order to prevent termination of simulation, either the time step had to be reduced or dx and dy had to be increased. Many different time steps were tried and failed. Coarsening the bathymetry grid solved the problem and significantly reduced the simulation run-time.
The next problem in the calibration was the abnormality in the velocity plots. Sudden jumps in velocity occur with immense oscillations during the transition from one month to the next, with values of subsequent months doubling with time [Fig 4.1a]. Different approaches were tried in bid to reduce the rise in velocity.
Fig 4.1a Velocity Jump between January to February
The initial effort was spent straightening the nooks and bends along the inflow tributaries to smoothen the flow and in attempt to ease the velocity. Speculating that the jumps between each month might be due to the large difference in monthly inflow data, values during the transition period between each month were smoothed out [Fig 4.1b]. Besides that, efforts were also spent splitting inflow at Tengah tributary into 3 smaller inflows, resulting in a total of 5 inflows.
Fig 4.1b Smoothed Inflow Data
These did not produce significant improvements and a further step was taken to trim the inflows [Fig 3c]. In doing so, the cross-sectional area is increased at the inflows, thus reducing the velocity for the same flow rate (V = Q/A). Velocity at the inflow would therefore be the maximum velocity due to the smallest cross-sectional area. Surface velocities should not exceed the velocity at the inflows as cross-sectional area increases along the reservoir.
Though slightly reduced in magnitude, the escalation of velocity continued with immense oscillations [Fig 4.1a]. In order to determine the driving force of the oscillations, two simulations with altered meteorological data were set up. The first simulation had no solar radiation, 100% relative humidity and a constant air temperature of 27.8 degC, ensuring that the only forcing is due to the effect of wind [Fig 4.1c]. In doing so, the only possible cause of evaporation is wind stress. The second simulation had no wind effects [Fig 4.1d], using real data for solar radiation, relative humidity and air temperature.
Fig 4.1c Pure Wind Stress Simulation
Fig 4.1d Simulation Without Wind Stress
Results show that the frequent oscillations were due to the effect of wind stress. Water level variations were significantly larger for the first simulation with wind stress compared to the second simulation. This shows that wind is the major driving force of hydrodynamics, affecting both surface water level and velocity.
Escalation in velocity was finally resolved after evaporation loss was catered for. ELCOM was initially modeled using mass balance where Σ inflows = outflow. This is impossible as evaporation was present. The effect of evaporation should cause the volume of water within the reservoir to decrease, in turn causing the water level to drop. Due to the mass balance specified for ELCOM, the velocities are escalating in order to maintain the amount of outflow with evaporation taking place at the same time. The outflow rate was reduced by an average of 15% to compensate for the evaporation loss and it significantly reduced the drop in water level to 0.15m for 6-months simulation [Fig 4.1f].
After calibrating ELCOM, the simulation spin-up time can be estimated to be 20 days. It can be observed in plots for tracer concentration [Fig 4.1e] and water level [Fig 4f] that the initial surge in values occurs over a period of approximately 20 days. After which, the trends are more gradual, likely due to evaporative losses and rather well-mixed conditions. Well-mixed condition is achieved when steady state is reached, where the upper and lower layers have duplicating trends.
Fig 4.1e Tracer Concentration against Time
Fig 4.1f Water Level against Time
The residence times of ELCOM can be obtained. Hydraulic detention time is the time taken for the water to reach the dam outflow from the inflow. It is estimated to be 15 months, given Ï„w = V/Q where V is the reservoir volume and Q is the sum of average monthly inflow from the three tributaries. This can be counter-checked using
Q = (annual rainfall x catchment area x runoff coefficient) / (365days x 86400s/d).
Residence time of tracer is the average length of time the tracer particle remains in the reservoir. It can be obtained from C = C∞ ( 1 - e-t/τ ). The gradient from the plot [Fig 4.1g] of ln(1 - C/ C∞ ) = (-1/τ)t using results from the 2-year simulation would produce the residence time of tracer. The residence time is 30.5 months, a multiple of the hydraulic residence time.
Fig 4.1g Tracer Residence Time
Upon establishing the tracer residence time, the following table can be determined:
Table 4.1a Tracer Concentration with Time
Tracer Concentration
Time
95%
91.3 months
50%
21.1 months
4.2 Results and Analysis
Results from ELCOM simulation were analysed, placing emphasis on temperature, tracer concentration, velocity and surface water level. In total, nine point profiles, three curtains, two sheets and a general 3D were produced.
The nine point profiles consists of the three inflow points, a point at the dam outflow, the deepest point within Kranji Reservoir and four PUB points (in blue on Fig 3c, labeled in the direction of flow). Three cross section curtains span along Kangkar (curtain 1, red-yellow), Tengah (curtain 2, yellow) and Pengsiang (curtain 3, green-yellow) to the dam outflow. The two sheets are the plan views of the top and bottom layer of the reservoir. This work shall concentrate on presenting the point profiles and curtain through result plots and movie snapshots.
Temperature plots for the three inflows show similar trends [Fig 4.2a]. The water column is rather well-mixed vertically as depth is shallow at about 1.5m. The top (blue) and bottom (red) layer generally overlap due to mixing. Slightly higher temperature in the surface layer is likely due to the effect of solar radiation. The simulated results show that the model can adequately simulate temperatures at surface, bottom and intermediate depths of the water column [Fig 4.2a(v)]. Temperature difference is more prominent in points of greater depth [Fig 4.2a(v), (ix)].
Fig 4.2a(i) Temperature at Kangkar Inflow
Fig 4.2a(ii) Temperature at Tengah Inflow
Fig 4.2a(iii) Temperature at Pengsiang Inflow
Fig 4.2a(iv) Temperature at Dam Outflow
Fig 4.2a(v) Temperature at Deepest Point
Fig 4.2a(vi) Temperature at PUB Point 1
Fig 4.2a(vii) Temperature at PUB Point 2
Fig 4.2a(viii) Temperature at PUB Point 3
Fig 4.2a(ix) Temperature at PUB Point 4
Temperature variations within Kranji Reservoir during different periods can be shown in movie snapshots of curtain 1 (Kangkar Infow). Fig 4.2b shows a typical warm day in the month of April in 2003 [Fig 4a]. The temperature of inflow is 29degC. It can be seen that the deepest point at 18m depth is cool compared to the surface temperature which is affected by solar radiation. As time passes, surface layers are warmed by the sun, reaching 33degC while the coolest bottom region remains at 28degC throughout the day.
Fig 4.2b Temperature Variation of Warm Day
Fig 4.2c shows the cooling down of the reservoir after it has sufficiently warmed up over a period of hot months[Fig 4a]. The bottom region has warmed up to 28.5degC but it still remains cool compared to the water mass of 29.5degC.
Fig 4.2c Cooling Down of Reservoir
Tracer (30units) was added in Tengah Inflow for the 2-year simulation.
The initial steep rise in tracer concentration for Tengah Inflow [Fig 4.3b(i)] is probably due to the initial spin-up of the model. Plots with less fluctuations is likely to be of deeper depth as condition is not as well-mixed. The transit time [Fig 4.3a] from PUB Point 1 (pink) to PUB Point 2 (blue) is approximately 6 days, 7 days from PUB Point 2 to PUB Point 3 (cyan) and 10 days from PUB Point 3 to PUB Point 4 (black).
Fig 4.3a Tracer Concentration of 4 PUB Points
Fig 4.3b(i) Tracer Concentration for Kangkar Inflow
Fig 4.3b(ii) Tracer Concentration for Tengah Inflow
Fig 4.3b(iii) Tracer Concentration for Pengsiang Inflow
Fig 4.3b(iv) Tracer Concentration for Dam Outflow
Fig 4.3b(v) Tracer Concentration for Deepest Point
Fig 4.3b(vi) Tracer Concentration for PUB Point 1 Fig 4.3b(vii) Tracer Concentration for PUB Point 2 Fig 4.3b(viii) Tracer Concentration for PUB Point 3 Fig 4.3b(ix) Tracer Concentration for PUB Point 4
The movement of tracer front through the reservoir can be presented by the movie snapshots of curtain 1 (Kangkar Inflow) [Fig 4.3c]. Movement of the tracer front is clearly visible as it gradually mixes and spreads throughout the water column. Tracer is mixing better within the initial inflow compared to the main water mass. After a month of simulation, the tracer concentration in the shallow region has reached an average of 10units while the deeper portion is still less than 5units.
Fig 4.3c Movement of Tracer Front
Surface water level fluctuations are within reasonable range of 0.3m after catering for evaporative losses [Fig 4.2d]. Mass balance Σ inflows = outflow + evaporation is achieved by reducing outflow rate by 15%. The plots for surface water level are almost identical for the nine different point locations within Kranji Reservoir.
Fig 4.2d Surface Water Level
Fig 4a 2003 Changi Meteorological Data Forcing Components
CHAPTER 5 - CONCLUSION
This study presents the assembly and calibration of ELCOM of Kranji Reservoir using meteorological data from Changi collected from January through December 2003. Grid resolution of 40 x 40m is sufficient to model the hydrodynamic processes within the reservoir reasonably well. The spin-up time for the model is estimated to be approximately 20days.
As shown in this work, wind effect is the primary driver of hydrodynamics. Regions of shallow water exhibit well-mixed conditions. The thermal balance is affected by air temperature and incident solar radiation. Tributary inflows and outflow have only local impact on the hydrodynamics.
Residence time can be obtained through the input of tracer. The tracer residence time is a multiple of the hydraulic retention time, which is a good indicator of the flushing time within the reservoir. Transit time from one point to another within the reservoir can also be estimated from the tracer concentrations.
The effect of rainfall has been neglected in this study. The simulation can be improved and modeled more realistically with the input of rainfall. This will affect the inflow within the overall mass balance and the evaporative loss.
REFERENCES
[1] Casulli and Cheng, 1992 - Semi-implicit finite difference methods for three-dimensional shallow water flow.
[2] Hodges, B.R, 2000 - Numerical techniques in CWR-ELCOM.
[3] Imberger and Hodges et al, 2003 - Modelling circulation in lakes: Spatial and temporal variations.
[4] Imberger et al, 2000 - Modelling basin-scale internal waves in a stratified lake.
[5] Leonard, B.P, 1991 - The ULTIMATE conservative difference scheme applied to unsteady one-dimensional advection.
[6] NTU Project Objective and Scope
[7] Simons, T.J, 1974 - Verification of numerical models of Lake Ontario.
[8] Spigel et al, 1980 - The classification of mixed-layer dynamics in lakes of small to medium size.
[9] The University of Western Australia, Centre for Water Research - www.cwr.uwa.edu.au
