## 1. Introduction

There are a number of directives to protect the environment, for example, in the EU see [1], [2]. Implementation of guidelines in accordance with these directives is based on environmental monitoring. Modelling based on existing data could look to reduce the necessary costs required for environmental monitoring in the future. This chapter was included in this book to reflect the adverse impacts of poor-treated and/or accidentally untreated water on nearby water bodies. Figure 1 shows an idealised sketch of mixing of pollutants in the river downstream of the outfall (0 m) of a wastewater treatment plant (WWTP).

In general four distinct zones exist downstream of the orifice. The first zone is a zone with background concentration. The second zone is a 3D mixing zone, the third zone is a 2D mixing zone and the fourth zone is a 1D mixing zone. In the initially three-dimensional mixing process (until the location where complete vertical mixing takes place) the drop-off of the maximum mass concentrationis relatively fast, *C*_{max}
*~ x*^{-1}, while it occurs more gradually in the vertically mixed,thus two-dimensional, mixing phase, *C*_{max}
*~ x*^{-1/2} [3]. Mixing in rivers can be described in terms of the near-field and far-field mixing zones. The near-field is defined as the area between the discharge structure and the location where the effluent is vertically mixed in the river. It encompasses buoyant jet mixing and boundary interactions. Far-field mixing occurs after the effluent plume is vertically mixed, and encompasses transverse and longitudinal mixing in the river due to buoyant spreading and passive diffusion [4]. 2D mixing zone modelling will be explained using a simple model in the MS Excel file named Czech immission test (CIT), which is prepared for determination of regulatory mixing zones located downstream of WWTPs in the Czech Republic. CIT is an MS Excel spreadsheet model of the solution to the two-dimensional advection dispersion equation. Therefore, it may only be applied at a far-field distance where the plume is completely or nearly-completely vertically mixed. An environmental impact assessment (EIA) is an assessment of the possible positive or negative impact that a proposed project (in our case the example WWTP) may have on the environment. The only effective way of monitoring changes in biological quality is to perform continuous biological monitoring with instruments suitable for early warning purposes, for example see [5]. Due to the random occurrence of accidental spills, the only way to quickly detect hazardous situations is to perform continuous monitoring of surface water quality. That’s the reason why the Daphnia Toximeter (a product of the German bbe Moldaenke company) has been installed at the monitoring station located on the borderline profile at the Odra River in the town of Bohumín. To monitor the biological impact of surface water quality on biota,this instrument uses young organisms of *Daphnia magna*. The organisms are exposed to the surface water in a flow-through chamber, to which the water from a monitored profile is pumped. New born daphnids aged less than 24 hours are placed into the chamber and stay there maximally for 7 days. They are then replaced by new ones. There are two reasons for this – older organisms are less sensitive than the young ones, and at the age of 9 to 11 days daphnids start to reproduce parthenogenetically, which leads to an undesirable change to the number of organisms in the chamber. This instrument can be checked through an Internet connection, thus making it possible to obtain online information about the biological impact of surface water quality in the monitored profile from anywhere with Internet access. Conversion of the daphnids’ behaviour to numerical parameters is computed. If there is an accidental leakage of a toxic substance above the monitoring station, the output is a toxic index alarm curve. The second part of this chapter describes an example of toxic index alarm curve analysis.

## 2. Problem statement

The advective-dispersive equation for solute movement through a river forms the basis of the mathematical algorithm used by the riverine component. The surface-water flow is assumed to be steady and uniform; the algorithms are developed for the limiting case of unidirectional advective transport with three-dimensional (longitudinal, lateral, and vertical) dispersion. The advective-dispersive equation for solute movement in a river can be described by the following expression

,where *C* is a dissolved instream contaminant concentration [kg m^{-3} or %]; *U* is the average instream flow velocity [m s^{-1}]; *Ex,Ey,Ez* are dispersion coefficients in the *x-, y*-, and *z*-directions, respectively [m^{2} s^{-1}]. If source terms *‘S’* and ‘*f*_{R}*’* are added as shown in the equation above, the so-called advection–diffusion reaction equation emerges. The additional terms represent [6]:

Discharges or ‘wasteloads’ (

*S*): these source terms are additional inflows of water or mass. As many source terms as required may be added to Equation (1). These could include small rivers, discharges of industries, sewage treatment plants, small wasteload outfalls and so on.Reaction terms or ‘processes’ (

*f*_{R}). Processes can be split into physical and other processes.

Examples of physical processes are:

settling of suspended particulate matter

water movement not affecting substances, like evaporation

volatilisation of the substance itself at the water surface.

Examples of other processes are:

biochemical conversions like ammonia and oxygen forming nitrite

growth of algae (primary production)

predation by other animals

chemical reactions on whether the flow is laminar orturbulent.

Dispersion is the scattering of particles or a cloud of contaminants by the combined effects of shear and transverse diffusion. Molecular diffusion is the scattering of particles by random molecular motions, which may be described by Fick’s law and the classical diffusion equation.Turbulent diffusion is the random scattering of particles by turbulent motion, considered roughly analogous to molecular diffusion, but with eddy diffusion coefficients (which are much larger than molecular diffusion coefficients). The diffusion coefficients would either be molecular or turbulent, depending on whether the flow is laminar or turbulent [7]. In natural rivers, a host of processes lead to a non-uniform velocity field, which allows mixing to occur much faster than by molecular diffusion alone [7]. Under the assumptions of negligible momentum and buoyancy, and for a discharge near the stream's free surface or near the bottom, complete vertical mixing is expected to occur at distance

,where *x*_{3D-max} is the distance to complete vertical mixing [m]; *U* is the mean velocity downstream of discharge outfall [m. s^{-1}]; *h* is the mean river depth downstream of discharge [m]; *E*_{z} is the vertical mixing coefficient [m^{2}.s^{-1}]; *β* is the vertical mixing coefficient constant [dimensionless]; *U*_{s} is the shear velocity [m.s^{-1}]; *g* is the gravity acceleration constant [m. s-^{2}] and *C*_{C} is the Chezy’s coefficient [m^{0,5}.s^{-1}].

The constant *β* can be derived from the velocity profile (for example, see [8]), and its value is recommended at approximately 0.07 ± 50%. The distance to complete vertical mixing will be somewhat reduced if the finite dimensions of the discharge opening are considered or if the source location is varied within the water column, e.g. located at mid-depth. Thus, the discharge design can play a certain role in this initial region [3].

Complete transverse mixing is expected to occur at a distance

,where *x*_{2D-max} is the distance to complete transverse mixing [m]; *U* is the mean velocity downstream of discharge outfall [m. s^{-1}]; *y*_{o} is the discharge distance from nearest shoreline [m]; *B* is the river mean width downstream of discharge [m]; *E*_{y} is the transversal mixing coefficient [m^{2}.s^{-1}]; *α* is the transversal mixing coefficient constant [dimensionless]; *U*_{s} is the shear velocity [m.s^{-1}] and *h* is the mean river depth downstream of discharge outfall [m].

The constant *α* is high in meandering and curved channels due to the presence of secondary currents. A value of 0.6 is recommended for most natural channels. Fischer [8] reports that this constant can range from 0.1 to 0.2 for straight artificial channels. Curves and sidewall irregularities increase the constant *α* such that in natural streams it is rarely less than 0.4. If the stream is slowly meandering and the sidewall irregularities are moderate, then the constant *α* is usually in the range of 0.4 to 0.8. Therefore, a value of 0.6 is usually recommended in natural channels. Uncertainty in this constant is usually at least ± 50 %.

In most practical problems we can start by assuming that the effluent is uniformly distributed over the vertical, or in the other words, we can analyse the two-dimensional spread from a uniform line source [8]. The longitudinal mixing term has very little influence on transverse mixing under the above condition and can be dropped. If the channel has the width *B* the effect of the boundaries, *y*=0 and *y*=*B*, can be accounted by the method of superposition described in [8]. Volume friction of effluent in a sample at some sampling point in the effluent plume is independent on the ambient background concentration. The above facts were used to construct the Czech Immission Test (CIT) model. This CIT model is based on the following equation

,where *VF*_{e} is the volume friction of effluent in a sample at some sampling point in the effluent plume [dimensionless]; *Q*_{e} is the effluent discharge flow rate [m^{3}.s^{-1}]; *U* is the mean velocity downstream of discharge outfall[m. s^{-1}]; *B* is the river mean width downstream of discharge [m]; *h* is the mean river depth downstream of discharge outfall [m]; *x´* is defined by setting *x´= xE*_{y}*/UB*^{2} [dimensionless]; *y´* is defined as *y/B* [dimensionless] and *y*_{o}*´* is defined as *y*_{o}*/B* [dimensionless].

Dilution factor *DF* is reciprocal of the *VF*_{e}.

The concentration of a pollutant of concern (or a dye tracer) at any sampling point is given by

,where *C*_{e} is the solute concentration in the effluent [kg.m^{-3} or %]; *C*_{r} is the background concentration in the receiving water [kg.m^{-3} or %]; *DF* is the dilution factor [dimensionless]; *VF*_{e} is the volume friction of effluent in a sample at some sampling point in the effluent plume; *k* is the first order reaction coefficient [s^{-1}] and *t* is the time [s].

The closeness of the approach of model values to measured values is given by:

,where *C*_{f} is the field value; *C*_{m} is the model value and *C*_{avgf} is the mean of the field values. The perfect fit is indicated by

The unsteady solutions of the 1D advective dispersion equation (ADE) may be obtained using methods of Fourier transform, of Laplace transform and the Fourier method. Using the methods of calculus, analytical solutions are developed that provide the predicted solute concentration as a function of time and space. Analytical solutions are derived for conservative substance, constant velocity, constant discharge, constant cross-sectional area and constant dispersion coefficient. Point sources such as accidental spills may be viewed as instantaneous sources [10]. For a spill, the solution to ADE can be obtained using the method of Fourier transform. For a continuous rectangular input the solution can be obtained using the method of Laplace transform described by [8]. This solution is useable only for the duration of the continuous input. The principle of superposition may be used to develop the solution for all time periods after the termination of the continuous input [10]. Results of Fourier transform have greater error then results from Laplace transform [11]. The analytical solution of Laplace transform gives results comparable to the more time and calculation consuming Fourier method [11].

Sometimes the concentrations in the 1D zone are too low. The International System of Units (SI) doesn`t include ppb as a unit. Therefore it is convenient to compute with ratio *R*_{C}:

,where *C* is the concentration at the station at time t [kg.m^{-3}] and *C*_{IN} is the input solute concentration after mixing over the cross-section of the stream computed specially for every monitoring station (effect of dead zones is included) [kg.m^{-3}]. Every sampling point has its own *C*_{IN} computed using Equation (8):

where *m*_{s} is the mass of substance present in the accidental leakage cloud as it passed the sampling point [kg]; *A*_{m} is the model cross-section area [m^{2}]; *U*_{r} is the retarded velocity (dead zones included) [m.s^{-1}];*τ* is the duration of the continuous input [s]; *x* is the distance between source and monitoring station [m] and *t*_{p} is the time to pick value [s].

(9) |

,where *x* is the distance between source and monitoring station [m]; *t* is the time [s]; *U*_{r} is the retarded velocity [m.s^{-1}]; *τ* is the duration of the continuous input [s]; *E*_{x} is the longitudinal dispersion coefficient [m^{2}.s^{-1}]; *erfc* is the complimentary error function and *K* is coefficient [dimensionless]. It was convenient to define dimensionless coefficient *K* for computation as follows: if t< τ then K = 0; if t> τ then *K* = 1. Numerous studies on longitudinal dispersion have been conducted over the past few decades. Fick's second law predicts how diffusion causes the concentration to change with time. In actuality, immobile-flow zones (dead zones) may invalidate Fick’s law. Chatwin [12] developed a method for determining longitudinal dispersion intended to address the problem of non-Fickian behaviour. Technically, the Chatwin’s method is only really valid for impulse releases, but it does provide a reasonable approximation for longitudinal dispersion for pulse and continuous releases [13]. Chatwin’s values *b* can be computed by equation

,where *b* is the Chatwin’s value [s^{0.5}]; *t*_{p} is the time to pick value [s]; *C*_{p} is the pick value [% or kg.m^{-3} or any other unit]; *t* is time [s] and *C* is value at value [% or kg.m^{-3} or any other unit].

Longitudinal dispersion is given by Equation (11)

where *b** is the vertical axis (at *t* = 0 s) intercept of the straight-line fit to the early-time values *b* [s^{0.5}] and *x* is distance [m].

## 3. Some existing models used in water modelling and model choosing

There are some predictive models for examining the mixing from point sources and showing compliance with EQS-values:

General water quality models may be required in more complex situations. Different methods for the far-field modeling exist, ranging from water quality models in estuary- type flows (e.g. model QUAL-2 of the U.S. EPA), to Eulerian coastal circulation and transport models (e.g. Delft3d of Delft Hydraulics) to Lagrangian particle tracking models. For rivers, the model AVG of the German ATV-DVWK, the model QUAL-2 of the U.S. EPA, and the model RWQM1 of the International Water Association are all examples of general water quality models. Such models also form the basis of management procedures for attaining a good quality status in the case of multiple sources, i.e. by following the principle of a distributed waste load allocation for individual water users [3]. Danish Hydraulic Institute MIKE Software is the result of years of experience and dedicated development. DHI Software models the world of water - from mountain streams to the ocean and from drinking water to sewage [14]. MIKE 11 is synonymous with top quality river modelling covering more application areas than any other river modelling package.

Choosing of an initial dilution model - five models are described in [15]:

Three are theoretical (UM, UDKHDEN, and VSW), and two are empirical (RSB and CORMIX). UM is the current version of the earlier models UOUTPLM and UMERGE. It acts as a two-dimensional model for single ports, though a pseudo-three-dimensional version is employed when there is a multiport diffuser with potential merging. It uses the 3/2 power profile to calculate the ratio and determine the centerline concentration as a function of the top hat concentration that it predicts. The ratio changes continuously with each integration step along the trajectory. Merging is simulated with the reflection technique. The CORMIX model has three modules: cormix1 for submerged single-point discharges, cormix2 for submerged multi-port diffuser discharges, and cormix3 for buoyant surface discharges.

Choosing a farfield model described in [15]:

There are two farfield models which are presently recommended for use. They are code named FARFIELD and RIVPLUM5. The appropriate farfield model to use in a particular mixing zone analysis depends on the combination of conditions involved:

The receiving water is sufficiently deep such that a plume will form and pass through the initial dilution phase without "Froude number less than 1", "overlap", or "boundary constraint" problems. Use FARFIELD as the algorithm (i.e., the version in 3PLUMES interface).

The receiving water is shallow and unidirectional; the effluent is thoroughly mixed surface to depth (i.e., no defined plume); and the discharge is a single port or short diffuser. Use RIVPLUM5.

There is/are bank constraint(s). Use RIVPLUM5, provided the conditions in 2. Above are also met.

Other shallow receiving waters (with no bank constraints) which occur with all other combinations of effluent plumes and discharger configurations. Use FARFIELD as a stand-alone model. A three-dimensional advective dispersion equation may also be appropriate.

The QTRACER2 (program for tracer-breakthrough curve analysis for tracer tests in karstic aquifers and other hydrologic systems)[13] is fast and easy method for evaluating tracer-breakthrough curves (BTCs) generated from tracing studies conducted in hydrologic systems. It has been reviewed in accordance with U.S. Environmental Protection Agency policy and approved for publication. Results may then be applied in solute-transport modeling and risk assessment studies.

## 4. Some details and results of 2D modelling

First, the model must be calibrated; that is, its parameters must be adjusted to match the behaviour of the prototype. Second the model must be validated. This means that a calibrated model must be compared to data not used in the calibration to determine whether the model is applicable to cases outside the calibration data set [7]. The CIT model was tested against data from a dye tracer study of the City of Arlington WWTP discharge to the Stillaguamish River, described in reference [16 ]. This study was performed on 22 August, 2006 and documented in a Mixing Zone Study report (CEG, November 2006, revised May 2007). The field study included injection of Rhodamine WT dye into the WWTP effluent at a known concentration; collection of bottled fluorescence samples from within the effluent plume; and measurement of river bathymetry, width, and current velocity. At seasonal low flow conditions observed during the dye study, the river was approximately 121 feet (36.9 m) wide with an average depth of 4 feet(1.22 m). Average current speeds, measured with a Swoffer meter, were 1.5 feet per second (0.46 m.s^{-1}). The river channel is relatively straight and uniform downstream of the outfall, and river cross-section bathymetry is similar at other locations up to 500 feet (152.3 m) downstream of the outfall. The outfall consists of a single port discharge (12-inch-diameter) discharging horizontally at the river bottom. The outfall discharge is located approximately 52 feet (15.86 m) from the left (south) bank at an invert depth of 4.61 feet (1.406 m) during low flow conditions. Appendix A contains plan and profile record drawings of the outfall. Effluent discharge flow through the outfall was 2.2 million gallons per day (8,318.4 m^{3}.d^{-1}) during the study. Manning’s roughness coefficient value (0.025) in the study report was based upon the average rock diameter observed at the site. Water column average dilution factors at plume centerline are summarized in Table 1. Except for the 50-foot distance, centerline profiles were measured over two time periods to better represent the time varying nature of the plume. The plume was observed to rise from the river bottom immediately following discharge to approximately river mid-depth, and was relatively unsteady with a billowing nature (wandering back and forth across river within a prescribed area). Between 100 and 304 feet from the outfall, complete vertical mixing of the plume was visually observed to occur, and the billowing nature of the effluent plume was less apparent. These observations are confirmed in the effluent volume fraction profiles at the mixing zone boundary (304 feet downstream), where both time period results are nearly indistinguishable, and effluent concentrations are nearly uniform from the top to bottom of the water column. Calibration of the RIVPLUM 5 model to the tracer study results produced a transverse mixing coefficient constant equal to 0.4.

Calibration of the CIT model to the tracer study results produced a transverse mixing coefficient constant equal to 0.408 (see Figure 2). The closeness of approach of model CIT values to field values was computed using Equation (6). Resulting value of CIT model was about the same as resulting value of RIVPLUM 5 model and within the range of experiments reported by [8], from which CIT and RIVPLUM5 were developed.

Although physical processes play a large role in determining the fate of solutes, chemical and biological processes may be equally important. We would like to describe results of modelling for non-conservative substances now. The next fictive exapmle shows possibility of investigating the impacts of different effluent quality using the CIT developer model to show its merits with different treatment efficiencies. We are interested about biological oxygen demand (BOD) in this case. We assume the calibrated and validated value of transverse mixing coefficient constant is equal to 0.408. Secondly we assume calibrated and validated value of the temperature dependent first order rate coefficient for BOD equal to 2.31.10^{-6} s^{-1} at a temperature 20°C. The CIT model computes concentration of BOD at point of interest using Equation (5). Resulting concentration of BOD at a downstream distance 200m can be found in Table 2.

Sudden changes usually occur in the regular channel of the natural river. These singularities can be natural, such as changes in roughness, meanders etc., or artificial, such as bridge piers, groynees etc. Burdych´s suggested method uses spare length *x*_{s} in these cases [17]. Spare length is a distance between the real location of discharge outfall and its fictive effective position. It is explained on the next fictive example. Effluent from the WWTP is discharged to the river at river kilometre 49.62. The river is approximately 50 m wide with an average depth of 3.32 m. The mean velocity of flow is 0.223 m.s^{-1}. The river is relatively straight and uniform downstream of the outfall. The outfall discharge is located at the bank of the river. Effluent discharge flow through the outfall is 0.455 m^{3}.s^{-1}. Field values of conductivity are listed in Table 3 and can be seen in Figure 3.

For each distance listed in Table 3 the value *F* was computed. We can distinguish three characteristic reaches. The average transverse mixing coefficient *E*_{y} for each reach can be obtained by fitting a straight line through values *F* plotted against the distance. *E*_{y} is the factor which determines the slope of regression line (see Figure 4). The intercept of this regression line can be expressed as *E*_{y}**x*_{s} [18].

The transverse mixing coefficient constant *α* can be computed using Equation (12)

,where *α* is the transversal mixing coefficient constant [dimensionless]; *E*_{y} is the transversal mixing coefficient [m^{2}.s^{-1}]; *C*_{C} is the Chezy’s coefficient [m^{0.5}.s^{-1}]; *h* is the mean river depth downstream discharge [m]; *U* is the mean velocity downstream of discharge outfall [m. s^{-1}]; *g* is the gravity acceleration constant [m. s^{-2}]; R is the hydraulic radius [m]; p is the Pavlovskij’s or Manning’s exponent and n is the roughness coefficient.

There are a few methods on, how to calculate Chezy’s coefficient, such as Manning’s method, Pavlivskij’s method etc. If we use Pavlivskij’s method *α* is equal to 0.199 in our case. If we use Manning’s method *α* is equal to 0.194 in our case. Another way to determine the average mixing coefficient constant is calibration using the CIT model (see Figure 5). The CIT model uses Manning’s method. Pavlivskij’s method is preferable for shallow rivers with great average rock diameter. Reach 1 needs negligible spare length, reach 2 needs the spare length equal to approximately 180 m and reach 3 needs the spare length equal to approximately 380 m.

## 5. An example of 1D modelling of accidental leakage

Thanks to the device for continuous monitoring of the biological water quality, which was installed at the Bohumín station located on the borderline profile at the Odra River, accidental leakage on the 9th of January 2012 was registered [19]. The first non-zero toxic index (*TOX*) was registered at 21:28:05. The maximal value of the toxic index was registered at 22:26:30. The difference between these two times is 3,505 s. The measured water discharge at the sampling station was 30.5 m^{3}.s^{-1}. There is usually a constant ratio between the time to pick of the value of the effluent cloud *t*_{p} and the time to beginning of the effluent cloud *t*_{b} at the far-field located sampling stations. This ratio is approximately between 1.1:1 and 1.3:1. You can find it in many studies based on tracer experiments, for example, see [20], [21]. In our case we find an optimal ratio of 1.17:1. The lag time of travel *t*_{L} is the time to the centroid of the toxic index in our case. One of the possible ways to compute *t*_{L} is shown in Table 4. The lag time of travel is particularly important because the centroid times of pollution concentration-time profiles at sites are often used to evaluate reach travel times (or velocities) for use in pollution incident and water quality models. Yet if the profiles so used are incomplete, the evaluation of the centroid by method of moments (the usual approach) is subject to error. Water and toxic pollution do not flow at one single advective velocity but experience a wide range of velocities, from rapid flow in the centres of large conduits to slow flow in adjacent voids. This variance of velocities is referred to as dispersion and is reflected in the width of the breakthrough curve.

In the second part of Table 4 the results for other characteristic times of effluent cloud passage at the sampling station are listed. The resulting lag time of travel (*t*_{L}=25,057 s) gives us distance between the source of accidental leakage and the sampling station of approximately 11,600 m. It is based on previous measurements of lag time of travel in the Odra River. An example of alarm curve is shown in Figure 6. The toxic index was calculated from the values of the parameters measured. The increasing part of toxic alarm curve has four characteristic values. The toxic index alarm curve is compared with the shape of the model R_{C} curve at Bohumín station on 9 January 2012. These solutions may be applied to our problem to determine the effects of advection and dispersion on solute concentration. Note that these equations do not consider any reactions the substance may undergo after entering the stream. As such, we are considering the worst case scenario, in which the toxic substance is transported conservatively and not allowed to degrade. This solution therefore provides an upper bound on the solute concentrations that are likely to be realized within the stream.

Chatwin’s values *b* were calculated using Equation (10). Results are listed in Table 5.

After fitting a straight line through the Chatwin’s values we obtained value *b** equal to 862.23 s^{0.5} (see Figure 7). The Equation (11) allows for solution of the longitudinal dispersion coefficient, provided that the plot against early-time data reasonably falls as a straight line. The late-time data will depart from the straight line due to non-Fickian dispersion characteristics (e.g., dead zones). For the distance of 11,600 m we calculated a longitudinal dispersion coefficient 45.25 m^{2}.s^{-1} using Equation (11). Ratio R_{C} was computed using Equation (9). Duration of accidental leakage *τ* was estimated so, that the shape of the *R*_{C} curve in Figure 6 is about the same as the shape of the toxic index curve. The potential source of accidental leakage could be at the bank of the Odra River approximately 11,600 m upstream of the sampling station, or at the bank of any tributary within this distance.

## 6. Conclusion

Water quality modelling can help us in environmental impact assessment. One of many 2D water quality models is the CIT model, which is described in this chapter. Every water quality model should be calibrated and validated for the conservative substance at first. Afterwards, this model may be calibrated and validated for non-conservative substance. The CIT model was successfully tested against data from the dye tracer study of the City Arlington WWTP discharge to the Stillaguamish River. CIT model may only be used at a far-field distance where the plume is completely or nearly-completely vertically mixed. Near-field mixing driven by jet velocity can be included to computing using spare length *x*_{s}. The buoyancy of discharge is neglected. Therefore the model conservatively underestimates mixing that occurs in the near-field if the plume is vertically mixed. The method described in the second part of this chapter could help us to find the source of accidental leakage, which is recorded at the sample station. The reason of this was to achieve about the same shape of the model R_{C} curve as the shape of toxic index alarm curve. Analytical solution of advection-dispersion equation was used. If numerical solution of advection-dispersion equation is used, then effect of numerical dispersion should be excluded. Described model might be used for early forecast of concentrations of conservative harmful substances in Odra River at Bohumín boundary profile.