HYDRODYNAMIC AND COHESIVE SEDIMENT TRANSPORT MODELING IN CHILIKA LAGOON

Chilika lagoon, one of the largest brackish water lagoons in Asia located along the east coast of India. The rivers draining into the lagoon carry about 13 million tonnes of sediments annually. Because of the cohesiveness properties of the fine sediments, nutrients, heavy metals and other polluted substances tend to bind to the sediment’s surface. Consequently, pollutants can be concentrated in the inlets/estuaries, thus being of great environmental interest. In addition, the mudflats occurring are important biotopes for a large number of microand macro-faunal species and act as feeding places for a number of birds. To understand the cohesive sediment dynamics, a numerical model, MIKE 21 Mud Transport (MT) coupled with hydrodynamic (HD) was used. The model simulated the relative bed level height and suspended sediment concentrations. The sediment interchange and accumulation between each sectors and Bay of Bengal were evaluated. The suspended sediment concentration is high in the north-east portion of the lagoon while medium and low suspended loads are observed in the eastern and western portion of the lagoon. Bed thickness is very high in the north-western corner of the lagoon covered with Phragmites Karka which facilitate sediment trap. Total bed thickness change is very much pronounced in the northern sector which receives most of the sediments from the Mahanadi river systems as well along the periphery of the lagoon due to drainage. The eastern lagoon shows a net deposition accumulated fraction (5-15 kg/m2) and hence gives enough indication of the sedimentation processes in the lagoon. Further, the results also warrant immediate attention to check and monitor suspended sediment concentration to find out the net deposition trend in the lagoon environment in order to take decisions in minimizing the sediment load.


INTRODUCTION
A coastal lagoon is a distinct dynamic environment where interplay of different energy forces from land-sea-atmosphere operates in a shallow body of water which is partly enclosed by a barrier and which has restricted or ephemeral communication with the sea through one or more inlets.The processes of flocculation, settling and scour lag, and the asymmetry of the tidal currents make cohesive sediment transport in the lagoon through the inlets difficult to forecast (Dyer, 1989;Teisson, 1991;van Leussen, 1994;Parker, 1997;Van der Lee, 2000).Because of the cohesive properties of the fine sediments, nutrients, heavy metals and other polluted substances tend to bind to the sediment's surface (Rae, 1997).Consequently, pollutants can be concentrated in the inlets/estuaries, thus being of great environmental interest.In addition, the mudflats occurring are important biotopes for a large number of microand macro-faunal species and act as feeding places for a number of birds (Eisma, 1998).This makes forecast of erosion, transport and deposition of cohesive sediment of great interest in all lagoon systems.Number of models is fairly well known in the hydraulic modeling community and has successful records of accomplishment on numerous applications throughout the world.Most of the present cohesive sediment transport models use an advection-dispersion equation to simulate the cohesive sediment transport in the water column (Teisson, 1997).The advection-dispersion equation requires current velocity components and water levels that are normally provided from a decoupled hydrodynamic model.The models have further incorporated a variety of equations that describe the cohesive sediment erosion, flocculation and deposition processes in different ways (Mehta et al., 1989;Teisson, 1997).Mahanty et al (2015) performed an analysis to identify the most appropriate hydrodynamic model to characterize the hydrodynamics and salinity of the Chilika lagoon.
Modelling the transport processes in a lagoon or in a marine environment (Platt et al., 1981;Petihakis, et al., 1999;Pinazo et al., 2004;Desmit et al., 2005) is of great interest because of its potential in preserving and predicting the ecology.Without a modelling framework and analysis, it is difficult to predict the alternatives available for decision making.Efficient modelling requires the knowledge of several related factors-dynamics / circulatory pattern (Blumberg, 1977;Rao, 1995;Rao et al., 1999) (transport by the moving waters), chemical (decay, reaction between substances), physical (transition between different states) and biological (migration of species), transport of nutrients.Dispersion of polluting substances (Fischer, 1976;Mazumdar and Das, 1992), sedimentation (Chandramohan et al., 1998;Mazumder and Ghoshal, 2002;Mazumder and Ghoshal, 2005;Mazumder and Ghoshal, 2006;Sinha et al., 2006), bottom erosion and species migration are all basically conditioned by the mixing and circulation of the water masses.The interaction of fresh water, saline water from the ocean and suspended sediments (Mazumder and Ghoshal;2002) adds to the mixing processes that are driven by the density differences.Considering the economical and Figure 1.Study Area environmental importance of cohesive sediment transport, accurate simulation of cohesive sediment transport processes in Chilika lagoon are necessary for environmental management, dredging studies, damming projects and even global sea level rise.

STUDY SITE
Chilika lagoon (19°28'N-19°54'N and 85°6'E-85°35'E) lies in the districts of Puri, Khurda and Ganjam in the state of Orissa, on the eastern coast of India (Figure 1).It is the largest brackish water lagoon in Asia with estuarine character.The lagoon is a highly productive ecosystem with rich fishery resources.The rich fishing grounds sustain the livelihood of more than 2,00,000 fisher folk who live in and around the Lagoon.Based on its rich biodiversity and socio-economic importance, Chilika Lagoon was designated as a Ramsar Site in 1981, especially as an important Water-fowl habitat.Hydro-logically, Chilika is influenced by three subsystems; the Mahanadi distributaries, the rivers/streams draining in to the lagoon from the Western catchment and the Bay of Bengal.Salinity is the most dominant factor determining the lagoon's ecology and the salinity dynamics are controlled jointly by the nature of the connection to the sea, associated tidal prism, and the volume and timing of freshwater inflow to the lagoon from the delta distributaries and western catchments.The vast amount of silt brought by Daya, Bharagavi and other streams contribute to the shallowness of the lagoon.Rao et al., (1986) made an initial geomorphic study of Chilika lagoon and adjoining area using remotely sensed data and showed that siltation process affecting the shrinkage of the lagoon area was significant.Water quality assessment of the Chilika lagoon using Indian Remote sensing Satellite (IRS) data was attempted by Sudhakar and Pal, (1993) and Pal and Mohanty, (2002) and the results showed significant temporal and spatial variability in the silt load.

Data
Bathymetry reflects the geometry of the region.The model bathymetries were prepared, based on the information from British Admiralty Sea Maps (extracted from a DHI C-MAP in digital form) and toposheets of Chilika region prepared by the Survey of India.The bathymetry map has been validated with the GPS observations collected during field survey in Chilika lagoon.All bathymetries used Chart Datum (CD), i.e. a datum set approximately equal to Lowest Astronomical Tide (LAT), and the projection used was lat/long.The mesh file for the Figure 2. The Flexible mesh structure of Chilika Lagoon Chilika lagoon (Figure 2) has been created considering the computational grid, water depths and boundary information.

Model Setup
Keeping in mind the need for a comprehensive sediment transport model for Chilika lagoon, an attempt have been made to setup a calibrated numerical hydrodynamic model (MIKE 21 HD) and further, to model cohesive sediment transport (MIKE 21 MT) in Chilika lagoon.The HD module calculates the resulting flow and distribution of salt and temperature subjected to a variety of forcing, sources and boundary conditions.The HD module solves the depth integrated incompressible Reynolds averaged Navier-Stokes equations (DHI, 2007a).The spatial discretisation of equations is performed using a cell-centered finite volume method.The spatial domain is discretised by the subdivision of the continuum into non-overlapping elements/cells.In the horizontal plane an unstructured grid is used comprising of triangle elements.2009), Lumborg and Windelin (2003) and Lumborg (2004).The module includes all important processes including dispersion of suspended sediment (Teisson et al., 1993), erosion from the sediment bed taking waves into account (Partheniades, 1965), settling of sediment using varying settling velocities (flocculation) (Krone, 1962;Burt, 1986;Pejrup, 1988a;Teeter, 1986), and deposition (Krone, 1962).MIKE 21 MT can take suspended transport of fine grained non-cohesive sediment into account.This is done by calculating an equilibrium concentration profile based on the sediment properties and the hydrodynamics.The bed is assumed to erode as flakes which mean that the distribution of fractions within the bed is also the distribution when eroded.This means that the erosion formula used in the MT section controls the maximum erosion of all fractions.After the flakes have been eroded it is assumed that they are destroyed or regrouped by turbulence.Since the sand fractions has no cohesive properties it will be freed by this and behave independently.The model does this by calculating the maximum possible equilibrium concentration for the given sand under the given hydrodynamic properties.If this is above the concentration of the sand fraction, the extra sand will be deposited so that the concentration is the equilibrium concentration.

Hydrodynamic Model
Hydrodynamics of the lagoon especially the circulation pattern is prerequisite to estimate the other processes in the lagoon.Circulation in the lagoon is governed by many forcing like bathymetry, wind stress, tides and freshwater influx from the rivers.In Chilika, wind driven circulation dominates the density driven circulation (Mohanty and Panda, 2009).Tidal influx also causes changes in the circulation in lagoons but its effect is limited to the region near the tidal inlet.Similarly, the effect of the freshwater influx contributes to the circulation in the main body of the lagoon apart from the tidal influence.For the present study, the hydrodynamic model was set up for two seasons (one month each in May and December 2009).
Calibrating a numerical model is an essential step in any modelling study.The observed field data was compared with the model results.In the hydrodynamic model (MIKE 21 HD) there are three calibration factors: the bed resistance coefficient, the eddy viscosity coefficient and the wind friction coefficient.Since there are not measurements for any of these factors, in Chilika lagoon, their values were defined during the calibration procedure.A Manning number of 32 m 1/3 s -1 was applied to the entire study area making the measured and the modelled water levels in good accordance throughout the study period.A good agreement found between simulated and observed water levels at Rambhartia (Figure 5a and 5b).

Sediment Transport
The cohesive sediment transport component of the model is generic and therefore requires several calibration factors to obtain a good description of erosion, transport, settling and deposition of the sediment.Erosion is initiated when a given bottom shear stress (the shear strength-τce) is exceeded.The description of the sediment bed was designed so that the shear strength is increasing with increasing depth as a dynamic bed description is required in long term simulations (Cancino and Neves, 1999).The shear strength τce-values in the interval 0.16-3 N m -2 have been used.The τce -values were varied down through the bottom layers so that the lowest values were used for the top layers.Further, at sites with more wave action and in the deeper parts where the bottom is less muddy higher τce -values were applied.Generally, the model sensitivity for τce-values is highest at inlets.The erosion formulation given by Parchure and Mehta (1985) further requires an erosion rate (E) and an erosion coefficient (α).These values determine the amount of sediment eroded per time provided that sce is exceeded and again the magnitude of these values have the highest impact on the intertidal flats.Sensitivity tests have shown that changes in E and a can change the net deposition profile significantly.The erosion coefficient (a) was used as a calibration factor and a value of 6.5 m N -0.5 was chosen for all bed layers.Transport of the suspended sediment is calculated using the advection-dispersion equation (Ekebjrg and Justesen, 1991;Teisson et al., 1993).The equation requires dispersion coefficients in the two horizontal directions.These coefficients are dependent on the mesh size and the calculation time step.In this study the dispersion coefficients were selected as 25 m s -1 .When the bed shear stress falls below a critical value (the critical bed shear stress for deposition-τcd) the suspended sediment will begin to deposit (Mehta and Partheniades, 1975).The settling velocity (ws) has earlier been shown to be of major importance when deposition is computed (van Leussen, 1994;Lumborg and Windelin, 2003).The settling velocity is dependent on several factors, the most important being the grain size, the suspended sediment concentration (Burt, 1986;Pejrup, 1988a), the turbulent shear stress in the water column (Manning and Dyer, 1999), and the biological activity in the system (Andersen and Pejrup, 2002).The relationship requires two site specific parameters and these have been obtained using settling tube measurements from the study area taken during December, 2009.The results yield the following algorithm: ws=3.96×10-6 ×SSC 1.19 , where ws is the settling velocity in ms -1 and SSC is the suspended sediment concentration in mg l -1 .The deposition is described using the deposition equation given by Krone (1962).The formulation is basically a settling flux giving the total settling as the product of the settling velocity (ws), the near bed sediment concentration (cb), and a probability factor (pd) which includes the critical shear stress for deposition (τcd).The formula performs well in the study area using low scd-values.After tests using different values, a differentiated solution was chosen with values in the range 0.05 to 0.3 N m -2 , using the lowest values in the Northen sector and higher values near outer channel areas of the lagoon.The order of magnitude is consistent with values found in the literature (Mehta and Partheniades, 1975;Cancino and Neves,1999;Whitehouse et al., 2000;Krishnappan and Marsalek, 2002).The concentrations at the open boundaries can be a crucial factor in a sediment transport model and even though the area of interest is located at some distance from the boundaries the concentrations here are selected carefully.The simulated distribution of Suspended Sediment concentration The suspended sediment concentration is high in the north-east portion of the lagoon while medium and low suspended loads are observed in the eastern and western portion of the lagoon (Figure 7a).Bed thickness is very high in the north-western corner of the lagoon (Figure 7b & Figure 7c) covered with Nala grass which facilitate sediment trap.Total Suspended Sediment concentration is below 0.75 kg/m3.In most parts of the lagoon (Figure 7d) only in few pockets the values exceed 0.5 kg/m3.Settling velocity of suspended sediments in kost parts of the lagoon are below 0.001 m/s (Figure 7e).Total bed thickness change (Figure 7f) is very much pronounced in the northern sector which receives most of the sediments from the Mahanadi river systems as well along the periphery of the lagoon due to drainage.Bed mass change are limited to 0-20 kg/m2 in most part of the lagoon (Figure 7g).The simulated spatial profile shows negligible erosion fraction (Figure 7h) while the net deposition fraction in some parts of the lagoon (Figure 8i) is quite noticeable (of the order of 0.008 to 0.04 kg/m2).The eastern lagoon shows a net deposition accumulated fraction (5-15 kg/m2) and hence gives enough indication of the sedimentation processes in the lagoon.Further, the results also warrant immediate attention to check the net deposition trend in the lagoon environment in order to conserve and preserve the lagoon for future generation.The model simulation appears to watch with the ground reality closely as far as sedimentation process is concerned.Therefore, it integrated with the required coefficients and data obtained through field experiments, the model can be used as a predictive tool in understanding the sediment transport and in preparing the sediment budget of the lagoon environment for its concentration and sustainable management.

CONCLUSIONS
Siltation is a major threat to the general ecology of the lagoon environment.Inflow of sediment to the tune of 13 million tonnes per annum via land drainage and their long resident time is the main cause of siltation.The problem of siltation is more acute in the Northern Sector and the Outer channel area.The Northern sector siltation is massive due to rapid growth of noxious weed and the sediment discharge through the riverine system.The siltation in the outer channel area affects the free movement of juvenile from the sea into the lagoon and vice-versa.As a result, loss of valuable species of prawns and mullets has been noticed during past few years (Samal, 2011).It has also been observed that the breeding and spawning ground of many important fishes, mollusks and crustaceans have been destroyed due to siltation.
Tide and Water Current at Rambhartia (Figure 1) were collected during May and September 2009 and used for model calibration.The tidal range varies between 32-94cm and 31-93cm during May and December 2009 respectively.The current near Rambharita ranges between 0.02 to 0.62 m/s during May 2009 and 0.001 to 0.51 m/s during December 2009.The fresh water discharge and sediment load into the lagoon have been inferred from Samal (2011).Higher discharge occurs during the months of June-September.Very less discharge is observed during nonmonsoon period.Figure 3 depicts the discharge into the lagoon during 2009.The maximum discharge through the Makara River (located at northern part of the lagoon catchment) during 2009 is estimated to be around 1729 million cumecs.Similarly, the minimum discharge of 103 million cumecs was observed in Kantabania River (located in western part of the lagoon catchment).Similarly, the sediment load into Chilika lagoon monitored by Samal (2011) depicted in figure 4. The Makara river contributes higher sediment load into the lagoon.

Figure 3 .
Figure 3. Water discharge into the lagoon during 2009

Figure 5 .
Figure 5. Validation of simulated surface water levels with observed surface water levels near new inlet mouth a) May, 2009 and b) December, 2009.The tidal amplitude is maximum at the inlet mouths of the lagoon which decreases gradually as we proceed inward from the inlet.It is seen that maximum water level variations recorded in Outer channel which is under the direct influence of tidal fluctuations.The amplitude of fluctuations is of very low order in the main body of the lagoon.Water currents in Chilika lagoon were simulated for flood and ebb conditions during May, 2009 and December, 2009.Figure 6a & 6b, depicts the simulated surface water currents in Chilika lagoon for flood and ebb conditions respectively during summer (May 2009).The maximum current velocity simulated near the new inlet mouth (0.18-0.20 m/s) and it is much lower in the other parts of the lagoon (0-0.02m/s).However, the water moves towards southern sector through the central sector with a speed of 0.02-0.04m/s (Figure 6a).Similar observation during ebb period (Figure 6b) shows that the Muggarmukha region which connects the outer channel with main body of the lagoon experiences highest water current of the order of 0.06-0.2m/s.The water current is very less (<0.02m/s) in northern sector and some parts of the southern sector.Major part of the central sector recorded relatively higher water current of order 0.02-0.04m/s.Similarly, Figure 6c and 6d represents the simulated water current for flood and ebb conditions during winter (December, 2009).The maximum flood condition water current is found between Muggarmukha to Satapada channel areas (0.18-0.28m/s) and gradually it decreases towards main body of the lagoon.The northern sector experiences the lowest water current (0-0.025m/s).The flood water enters into the lagoon through the Muggarmukha channel and diverges into two streams one is towards southern sector through central sector and other one towards northern sector.Further the northern sector stream diverges into two parts and forms eddy like circulation pattern in western and eastern parts of northern sector.Also the intrusion of rivers water from Nuna, Daya and Bhargabi rivers in north eastern side of northern sector plays an important role in circulation of this part as the directions of circulation are indicating.Similarly, during ebb period, the Muggarmukha to Satapada upto inlet mouth experience highest water current (0.125-0.25 m/s).The water recedes with a speed of 0.25 to 0.1 m/s from southern sector towards channel areas through the central sector.The rest of the lagoon including whole of the northern sector experiences very low water current of the order 0 to 0.025 m/s.The depth gradient results the stronger ebb current than of flood current.Experiments show that the tides are mostly effective around the channel area.The new inlet mouth opening has helped in increasing the tidal influx and hence the salinity, its influence is still not felt far interior of the lagoon due to the constriction of flow area between the lagoon and the channel area near Muggarmukha.