Analytical model of heavy-particle pollutant transport in estuary

The aspects of heavy-particle diffusion in the tropical estuary is investigated. The two dimensional advection-diffusion equations with tidal and settling velocity were solved with steady state approximation. An analytic solution was obtained in term of Cosines Integral function. The results showed that the profile of concentration follow a tidal cycle with the concentration decrease when the tidal speed is low and increase as the tidal speed is up. On the other hand, at the low speed of tidal current, the pollutants will likely to stay in place.


INTRODUCTION
Most coastal areas around the world are characterized by the presence of estuaries and embayment in which the water is predominantly driven by tides.The tide, especially tidal current is an important factor in dispersion processes of the passive scalar (pollutant) in the estuary.In general, ocean currents generated by tidal motion do not have a significant influence in the current circulation off the coast.In the area near the coast a tidal current has a significant effect even though we only discuss behavior to pollutants dispersion in the steady state condition.In the area near the coastline or estuary, the tidal current strengthened due to the shoaling effect.This flows generate the mixing processes and stratification in shallow waters which is closely related to the dispersion of pollutants.In short, the tide current is very important to support the pollutant dispersion in the estuary and the shallow coast.
A shallow tidal coast in which the ocean dynamics were dominated by tide is one of the most complex coastal seas.For example, the tidal current can generate the high level turbidity (Hoitink, 2003).The asymmetry of tidal currents is also a very important factor in the transport and accumulation of sediments (Hoitink, 2003;Wolanski et al., 1996).The works show that asymmetry in the flow acceleration will increase the rate of suspended load transport for fine sediments.The observation indicates that tide-dominated estuaries usually characterized by a convergent shape.
The complex processes of transport of pollutants in coastal areas is not only determined by a tidal hydrodynamics, but also by the variable propensity of many such contaminants (Eric et al., 2013).The example of this problem is that we are dealing with heavy-particle pollutants such as a heavy metal (Cd, Cu, Hg, Pb, Zn).The heavy metal measurement performed on the Scheldt River basin show that a longitudinally dissolved Cd and Cu has a convex shape instead Gaussian as most pollutants with fine particle (Gao et al., 2013).The studies suggest there are three things that need attention in the study of heavy metal dispersion.They are, inhomogeneity of point source for different metals, difficulty in quantifying diffusive input and complication of the biogeochemical behavior.Some effort to understand the process is done by using mathematical models such as 1D dispersion model, WASP water quality model, Environmental Fluid Dynamics Code (EFDC), SiAM-3D etc (Eric et al., 2013;De Smedt et al., 2006, Lu et al., 2014;Thouvenin et al., 2007).
In line with numerical model, the analytical model is also worth studying due to the simple nature.The advantage of an analytical model is simplicity so that the solution can be easily applied to field observations.The model usually only consider the importance or dominant factor in the phenomena such as estuaries dominated by dispersion, dominated by advection, or having dispersion and advection of similar magnitude (Fernando et al., 2012;Savenije 2005).For example the analytic solution of a simple depth-averaged 1D model has been used to investigate the sediments transport measurement with only consideration of the advection due to tidal motion, local resuspensions and deposition (Yu et al., 2012).
The paper investigates the analytical modeling of the heavy-metal dispersion in an estuary.In this proposed model, the heavy metal is considered as heavy-particles where its gravitational effect significant enough and this must be taken into account in the model.This effect will be accommodated by the settling velocity ws.The settling velocity is the velocity at which pollutant falls through water and it is controlled by both the drag force and the gravitational force.The combined effect of gravitation and the tidal current is challenging to model in analytical ways.Although the analytic model can only be done in the simple case, it is still important because the physical interpretation is easy to be obtained.In this paper we will investigate how heavy-particle pollutants disperse in the estuaries by solving an advection-diffusion equation including settling velocity and tidal currents.The steady state of two dimensions of an advection-diffusion equation is solved by using the separation variable methods.The analytic solution with periodic advection velocity and polynomial diffusion coefficient were obtained in term of Cos Integral function.
The paper is organized as follows: First is a description of heavy-particle pollutant dispersion followed by that of the analytic solution based on the separation variable.The advection processes due to tidal dynamics were next derived; thereafter, the simulation result of pollutant dispersion due to variation of horizontal and vertical coefficient were presented respectively, followed by the Mubarak and Sulaiman 41 concluding remarks.

PROPOSED MODEL
Pollutant dispersal around the estuaries will depend on two factors, that is, advection processes due to hydrodynamics motion and diffusion factors due to turbulence motion.By assuming advection velocity flow faster than pollutant dispersion, pollutant dispersion can be viewed as a steady state.Therefore, the steady state of an advection-diffusion equation describing dispersion of heavy-particle in the estuary can be written as, (1) where Kx and Kz are diffusion coefficients in zonal and vertical direction respectively, U is the zonal velocity and ws are the particle settling velocity.The settling velocity is the key factor to study the suspended load transport where the response of the sediment is characterized by H/ws with H denoting the water depth.It is well known that increasing of flow rate can entrain the sediments near the bed and then moves vertically upward until a steady state is reached (Hoitink, 2003).Previous studies show that most of the solution to be found are for special cases of this equation (James, 2002;Stockie, 2012).The boundary value problem is depicted in Figure 1.In this paper a special case will be used to get an analytic solution.The boundary and initial condition for C(x, z) are written as, (2) For simplicity, first it was assumed that diffusion coefficients are constant, say K.It is not difficult to show that generalization can be inferred from the solution, if the diffusion coefficients are not constant.Further, let us assume that the depth of estuary is constant.This approximation is useful for the alluvial estuaries.
In the coastal ocean, there is a disparity between horizontal and vertical scales.It means that the horizontal and vertical diffusion are usually considered separately and have very different values (James, 2002).It can be inferred that, an assumption that the separation variable is a solution can be made.By using the approximation, suppose the solution has the form of, (3) The equation can be written in two ordinary differential equations as follow, Next we solve the second equations with wrote back into, In the first step, it must be shown that the differential equation is exact.The equation is exact if we choose m = 0.This means that the equation should be solved, An exact linear second-order ordinary differential equation (ODE) can be solved by reduction to a linear first-order ODE (Kreyzig, 1993).This equation can be integrated if the dispersion coefficient is prescribed.For the most case in the estuary, the vertical motion is dominated by mixing processes due to the meeting of fresh water from a river and salts water from the off-shore.The mixing also caused by the breaking of internal waves due to shoaling effect when the wave propagates into the coastline.The physical mechanism of the processes is called the Kelvin-Helmholtz instability where the detail study is out of this paper (Smyth and Moum, 2012;Thorpe, 2007).The mixing experiment in the estuary has been done to get the vertical dispersion coefficient (Thorpe, 2007).The result shows that in the most case the vertical dispersion coefficient has the power series form or logarithms profile especially for close to the seabed.Following this result, it was assumed that the vertical dispersion has the form of K(z) = az b with the solution being, This is obtained by using DSolve Mathematica Package.
Further, the equation for ϕ(x) should read, This is not an exact differential equation as such cannot be reduced to the first order differential equation.That equation can be written in the form of the following, There is no general solution of this equation that must take special form of Kx(x) and U(x).To the pollutant dispersion in an estuary, the dispersion coefficient generally can be approached by the form of a polynomial function which is Kx(x) = px q with p and q are constant.In the alluvial estuaries the longitudinal dispersion usually has the form of exponential (Vallino and Hopkinson, 1998).Interest in the dynamics of pollutants due to tidal current will result in the flow taking the form of function U(x) = An cos(nπx).The general solution is given by, where Г(ζ) is a Gamma function.For example in the case of Kx(x) = px the solution is, where Ci(ζ)=γ+ln(ζ)+∫0 ζ (cos(η)-1)/ηdη is a special function called Cos-Integral and γ = 0.5772156 is Euler's constant.Now the complete solution is given by, with the boundary condition given by Equation 4.

ADVECTION PROCESSES BY TIDAL CURRENT
An equation of motion that describes the dynamics of tidal elevation and tidal current is the shallow waters equation.This research is only interested in the movement of the dynamics of the zonal direction that is dominant in the estuaries and near the coast, the dynamics of tidal elevation and current satisfy, with η and U are long wave elevation and particle velocity respectively.H is the water depth.For the nonlinear term is not dominant, the perturbation methods can be applied.In this approximation, the tidal elevation η and velocity U were described as, where ε= H/L is a small parameter and L, the zonal length scale.
Note that in shallow water approximation H ≪ L was used.The physical meaning is that the higher order has the small influence of the smaller order; in otherwords the higher order act as a correction of the lower order.Substituting Equations 17 and 18 into Equations 15 and 16 yields, The first orders have the solution, Where O is the term of third order and higher, An and Bn are the wave amplitude respectively.
The dispersion relations must be fulfilled is, Further, the second order can be written as, This is the inhomogeneous wave equation where the solution is given by (Officer, 1978), With the same manner, we obtain U (2) as, This is a simple form of tidal elevation and current.In the simulation we use steady state tidal current as follow, where this current can be obtained from a tidal measurement by applying the dispersion relation Equation 22.

RESULTS AND DISCUSSION
Here an example is given of a simulation model that was previously developed.The study areas to be learnt includes Tanjung Butung and Sumatra Indonesia.The estuaria is characterized by a micro tidal and gently sloping bathymetry.The time series of tidal elevation of this region is depicted in Figure 2. The first step in calculating pollutant dispersion, is to specify the tidal decomposition.Tidal phenomenon and the dynamics of sea associated with it in this region and the Indonesian water are very complex as the result of numerous strait, the complex topography as well as many small islands.Diurnal tides are unusually strong and are dominant along some coastlines (Ray et al., 2005;Gordon 2005;Wyrtki, 1961).From the figure, it can be seen that the tidal type is the mixed predominantly diurnal form.The tidal composition can be obtained by applying the Fourier transform respect to time series data.This is depicted in Figure 3.The tidal compositions from high to low energy are 22.98, 24.65, 11.78 and 11.50 h respectively.Substituting these periods into Equation 23 with x = 0 yields the tidal currents that is depicted in Figure 3.The advection velocity U(x) is obtained by substitution of tidal frequencies into the dispersion relation (Equation 22).With this relation, the wavelengths associated with its periods are 324.98, 348.61, 166.59 and 162.63 m, respectively.The tidal elevation associated with the wavelength is given by, where we used the cosines function.Then the tidal current then can be obtained by substituting into Equation 27.Finally, the tidal current decomposition into its component is depicted in Figure 4.The settling velocity 0.2 to 0.5 ms -1 was used.The next step is determining the diffusion coefficient for .The dispersion coefficient is larger to the offshore.This result is also confirmed in the study of Christensen and Li (2014).Generally, the studies show that for oceans and lakes, the horizontal dispersion coefficient (Eddy diffusivities) κ x and κ y are in the range of 10 4 to 10 8 cm 2 s -1 , and the vertical diffusivity κ z is between 10 1 to 10 3 cm 2 s -1 (Christensen and Li, 2014).In this paper, many cases of the form of longitudinal and vertical dispersion coefficients have been described.In the first case, the constant form of longitudinal dispersion coefficient namely K x (x) = p was chosen.The vertical dispersion coefficient is chosen in the form of K z = κ.The solution is given by, By using the Taylor expansion of this yields, The result is depicted in Figures 5 and 6.This is the pollutant concentration at the depth of 5 m at one tidal circle with the tidal wavelength about 350 m.Generally, the pattern of pollutant concentration is similar to the tidal current where the high concentration increases at high tide and low tide.A variation of the concentration occurred at the time ahead of high tide.The concentration profile showed that in general, the concentration will decline with the depth.The varying concentration occurred in the surface until 20 m depth.This means that the upper part of estuaries is act as the net of pollutant traps.
The second case, we choose the linear function for the form of longitudinal dispersion coefficient namely K x (x) = px.The vertical dispersion coefficient is chosen in the form of K z = κ.The solution until second order accuracy is given by, The result is depicted in Figure 7.The figure showed the pollutant diffusion in a tidal circle at five meters in depth.
The profile of concentration in general is still follow a tidal cycle.The concentration decrease when the tidal speed is low and increase as the tidal speed up.The pollutant diffusion in the water body is depicted in Figure 8.The figure showed that the concentration of going down to the depth with a pattern follow a tidal cycle.At the time of low tide phase, pollutants will likely to stay in place.
In the other case, it is assumed that there was a nonlinear form of the longitudinal dispersion and linear for the vertical dispersion coefficient, that is, K(x) = px 2 and K(z) = κz.Then the solution until second order accuracy is given by,

Conclusion
An analytical model of heavy-metal dispersion induced tidal current has been investigated.By applying separation variable methods the analytic solution was obtained.The explicit solutions can only be obtained by selecting particular form of the horizontal coefficients dispersion.In the case of dispersion coefficient, the linear and quadratic solution is given in the Cos Integral and Sin Integral function.This research studied the diffusion of pollutants in a tidal circle.The profile of concentration in general follows a tidal cycle with the concentration decrease when the tidal speed is low and increase as the  tidal speed up.All cases show that at the low speed of tidal current, the pollutants will likely stay in place.The expansion of the model is necessary for the case to be more realistically indispensable as taking into account the effects of the geometry of estuary.Previous studies show that the geometry is one of the most important parameters in the analytical estuary model (Tseng, 2002;Nguyen et al., 2012).On the other hand, the analytical model of the pollutant transport due to the residual currents produced by mixing of fresh water and salt water is also interesting and this work is still in progress.

Figure 1 .
Figure 1.The boundary problem of heavy-particle dispersion in estuary.The horizontal axis is represented the distance from offshore to the river mouth.
√gH is called linear tidal speed.Hence the tidal current was obtained by completing Equation 19 and 20 as follow,

Figure 2 .Figure 3 .
Figure 2. Tidal elevation from Tandjung Butung, Sumatra with time interval in hours and elevation in meters.This is a tidal prediction with time period March, 18 2015 to May, 18 2015

Figure 4 .
Figure 4.Tidal current decomposition for one tidal period.

Figure 5 .
Figure 5.The normalized pollutant concentration (C/ C0) in one tidal circle at z = −5 m for the constant horizontal coefficient.

Figure 6 .
Figure 6.The two dimensional dispersion of normalized pollutant concentration (C/C0) in one tidal circle with the horizontal dispersion coefficient is a constant.

Figure 7 .
Figure 7.The normalized pollutant concentration (C/ C0) in one tidal circle at z = −5 m with the linear horizontal coefficient.

Figure 8 .
Figure 8.The two dimensional dispersion of normalized pollutant concentration (C/ C0) in one tidal circle with the linear horizontal dispersion coefficient.

Figure 9 .
Figure 9.The normalized pollutant concentration (C/ C0) in one tidal circle at z = −5 m for the nonlinear horizontal coefficient.

Figure 10 .
Figure 10.The two dimensional dispersion of normalized pollutant concentration (C/ C0) in one tidal circle with the nonlinear horizontal dispersion coefficient.