**Modelling of Hydrobiological Processes in Coastal Waters**

Dimitrios P. Patoucheas and Yiannis G. Savvidis

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/59498

#### **1. Introduction**

It is well known that the dispersion of organic and inorganic matter in a coastal basin is closely related to the circulation of the seawater masses. The more detailed, accurate and reliable the knowledge of the hydrodynamic circulation is, the more successful the initial tracing of matter distribution in a coastal basin is expected to be.

In general, marine hydrodynamics is related to many different physical, chemical, geological and biological processes, which often constitute important environmental problems that need special care and detailed investigation. Among the most frequent phenomena of the afore‐ mentioned environmental issues, observed mainly in coastal basins, are the algal blooms, which are strongly related to optical, aesthetic and other important disturbances, such as eutrophication processes. The appearance of such phenomena is becoming more and more important, because algal blooms may contain potential toxic populations, called harmful algal blooms (HAB). As such episodes may correlated with eddies development [1], their study may obviously furnish a potentially effective prediction tool. The next important step in this rather preliminary process of tracing the potential locations of harmful algal blooms is the develop‐ ment and application of a transport model, which computes the matter transfer in space and time. A further step is the integration of a biological model into a transport model. It is expected that mathematical models that include physicochemical and biological processes will be able to describe the marine hydrodynamics and the relevant matter transfer issues, such as episodes of the spread of phytoplankton cells, contributing to the better understanding and more effective investigation and management of similar phenomena. Such models have been applied and discussed in literature [2-6].

The models to be used depend on the available data, the specific scientific questions to be answered and the particular characteristics of the ecosystem. The mathematical models used

© 2015 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

for hydrobiological processes in coastal waters can be focused on diagnostic matters, like the trophic state and the water quality of a specific coastal zone, or on prognostic matters, like the dispersion of HAB in a coastal basin.

Consequently, the use of mathematical simulation is of great significance for the diagnosis as well as the prognosis and prevention of hazardous situations that may threaten the marine environment. From this point of view, the description and application of three different types of models, presented in the following sections, constitutes the objective of the present chapter.

#### **2. Methods and case studies**

#### **2.1. Trix — Estimating the trophic status of a coastal ecosystem**

#### *2.1.1. Description of the method*

The trophic index (TRIX) proposed by Volleinweider et al. [7] is estimated as a linear combi‐ nation of the logarithms of four variables, referring to primary production [chlorophyll-a (Chla), dissolved inorganic nitrogen (DIN), total phosphorus (TP) and the absolute percentage deviation from oxygen saturation (aD%O)], in the form:

$$TRIX = \frac{\log\left(Chl - a\*aD\%O\*DIN\*TP\right) + k}{m} \tag{1}$$

Using data from the Northern Adriatic Sea, Volleinweider et al. [7] suggested the following values for the parameters *k* and *m*: *k*=1.5 (to fix the differences between the upper and lower range of the variables to 3 log units) and *m*=1.2 (to fix the scale range from 0 to 10). Low TRIX values (2-4) indicate poor productivity in waters, and thus imply "high quality" waters, while TRIX values from 4 to 5 indicate "good quality", 5 to 6, "mediocre quality" and 6 to 8, "bad quality" waters [8].

In this form, TRIX has been applied for the classification of coastal waters in the Mediterranean Sea, the Marmara Sea, the Northern European seas, the Black Sea, the Caspian Sea, the Baltic coastal waters, the seas of Southeast Mexico, and the Persian Gulf [9-17]. However, when the concentrations of the variables appear with different upper and lower ranges in respect to those of the Northern Adriatic Sea, the TRIX index needs modification by re-estimating the *k* and *m* parameters according to the new upper and lower limits [7, 18].

An unscaled trophic index, UNTRIX, was introduced by Pettine et al. [10] following the European Water Framework Directive (WFD, 2000/60/EC):

$$\text{LINTRIX} = \log\left(\text{Chl} - a \ast a \text{D\%} \prime \text{O} \ast \text{DIN} \ast \text{TP}\right) \tag{2}$$

Two classification procedures for water were proposed. The first is based on a simple com‐ parison between box and whisker plots for UNTRIX data from: 1) the investigation site and 2) the respective data from the reference area. The second procedure is based on the TQRTRIX trophic index, defined according to equation (3):

$$\text{TQR}\_{\text{TRX}} = \overbrace{^{50^{\text{th}} \text{UNTRIX}} \text{X}\_{\text{reference}}}^{\text{( $\text{t}$ )} \text{UNTTRIX}\_{\text{site}}} \tag{3}$$

where 50th is the median and 75th *percentile* is a value indicating that 75 % of the values of the variable fall below that value.

#### *2.1.2. Case study*

for hydrobiological processes in coastal waters can be focused on diagnostic matters, like the trophic state and the water quality of a specific coastal zone, or on prognostic matters, like the

Consequently, the use of mathematical simulation is of great significance for the diagnosis as well as the prognosis and prevention of hazardous situations that may threaten the marine environment. From this point of view, the description and application of three different types of models, presented in the following sections, constitutes the objective of the present chapter.

The trophic index (TRIX) proposed by Volleinweider et al. [7] is estimated as a linear combi‐ nation of the logarithms of four variables, referring to primary production [chlorophyll-a (Chla), dissolved inorganic nitrogen (DIN), total phosphorus (TP) and the absolute percentage

*m*

Using data from the Northern Adriatic Sea, Volleinweider et al. [7] suggested the following values for the parameters *k* and *m*: *k*=1.5 (to fix the differences between the upper and lower range of the variables to 3 log units) and *m*=1.2 (to fix the scale range from 0 to 10). Low TRIX values (2-4) indicate poor productivity in waters, and thus imply "high quality" waters, while TRIX values from 4 to 5 indicate "good quality", 5 to 6, "mediocre quality" and 6 to 8, "bad

In this form, TRIX has been applied for the classification of coastal waters in the Mediterranean Sea, the Marmara Sea, the Northern European seas, the Black Sea, the Caspian Sea, the Baltic coastal waters, the seas of Southeast Mexico, and the Persian Gulf [9-17]. However, when the concentrations of the variables appear with different upper and lower ranges in respect to those of the Northern Adriatic Sea, the TRIX index needs modification by re-estimating the *k*

An unscaled trophic index, UNTRIX, was introduced by Pettine et al. [10] following the

*UNTRIX = Chl a aD%O DIN TP* log( -\* \* \* ) (2)

and *m* parameters according to the new upper and lower limits [7, 18].

European Water Framework Directive (WFD, 2000/60/EC):


log(*Chl a aD%O DIN TP + k* ) *TRIX =*

dispersion of HAB in a coastal basin.

24 Hydrodynamics - Concepts and Experiments

**2. Methods and case studies**

*2.1.1. Description of the method*

quality" waters [8].

**2.1. Trix — Estimating the trophic status of a coastal ecosystem**

deviation from oxygen saturation (aD%O)], in the form:

The three abovementioned trophic indices (TRIX, UNTRIX and modified TRIX) were applied to quantify the water quality of the coastal waters in the Kalamitsi area, Ionian Sea, Greece (Fig. 1), before and after the operation of the Waste Water Treatment Plant (WWTP) [11]. As mentioned, this quantification allows a classification of the water quality.

Figure 1. The study area of Kalamitsi, Ionian Sea, Greece (Google Earth) **Figure 1.** The study area of Kalamitsi, Ionian Sea, Greece (Google Earth)

data set, following this procedure:

This procedure furnished equation (4):

method

deviation of the sample)

Data were collected seasonally, from three fixed stations, from March 2001 until January 2003 (before WWTP operation) and from August 2004 until July 2006 (after WWTP operation). Before the WWPT operation, the TRIX values showed a temporal variation in the range of 1.9-4.7 Data were collected seasonally, from three fixed stations, from March 2001 until January 2003 (before WWTP operation) and from August 2004 until July 2006 (after WWTP operation).

TRIX units, while no spatial variation (between stations) was observed. The mean TRIX value was 3.5. After the WWPT operation, the TRIX values showed temporal and spatial variations (2.2-4.8 TRIX units). The mean TRIX value was 3.4. As can be seen, no significant differences in TRIX values were detected before and after the WWPT operation (T-test, p=0.608). These values indicate "high" water quality, in contrast with the TRIX values obtained for the gulfs of the Aegean Sea, like the Thermaikos Gulf (5.0-6.0) and the Saronikos Gulf (3.7-6.2), with a mean of 5.3, indicating "mediocre" water quality [14]. As the concentrations of Chl-a and DIN from Kalamitsi displayed a profile much lower than those of the Northern Adriatic Sea [11, 7], the TRIX was re-estimated (KALTRIX) according to the Kalamitsi Before the WWPT operation, the TRIX values showed a temporal variation in the range of 1.9-4.7 TRIX units, while no spatial variation (between stations) was observed. The mean TRIX value was 3.5. After the WWPT operation, the TRIX values showed temporal and spatial variations (2.2-4.8 TRIX units). The mean TRIX value was 3.4. As can be seen, no significant differences in TRIX values were detected before and after the WWPT operation (T-test, p=0.608). These values indicate "high" water quality, in contrast with the TRIX values obtained for the gulfs of the Aegean Sea, like the Thermaikos Gulf (5.0-6.0) and the Saronikos Gulf (3.7-6.2), with a mean of 5.3, indicating "mediocre" water quality [14].

3



KALTRIX = Chl a aD%O DIN TP 2.5 log ∗ −∗ ∗ ∗ − ( ) (0.5) (4)

KALTRIX values displayed high temporal variability (0.8-9 units) before the WWPT operation, whereas both temporal and spatial variability occurred after the WWPT operation. In both cases, KALTRIX values were higher than those of TRIX, although the same data set was applied. The mean KALTRIX value was 5.3 before WWPT operation and 5.0 after it, with no statistically significant differences (T-test, p=0.735), indicating minor influence on water quality of WWTP operation. However, greater divergences of KALTRIX values were observed after WWTP operation. So, it could be assumed that KALTRIX is more sensitive, and is able to discriminate minor trends not easily



As the concentrations of Chl-a and DIN from Kalamitsi displayed a profile much lower than those of the Northern Adriatic Sea [11, 7], the TRIX was re-estimated (KALTRIX) according to the Kalamitsi data set, following this procedure:


This procedure furnished equation (4):

$$\text{KALTRIX} = 2.5 \ast \log \left( \text{Chl} - a \ast a \text{D\%} \text{O} \ast \text{DIN} \ast \text{TP} \right) - \begin{pmatrix} 0.5\\ \end{pmatrix} \tag{4}$$

KALTRIX values displayed high temporal variability (0.8-9 units) before the WWPT operation, whereas both temporal and spatial variability occurred after the WWPT operation. In both cases, KALTRIX values were higher than those of TRIX, although the same data set was applied. The mean KALTRIX value was 5.3 before WWPT operation and 5.0 after it, with no statistically significant differences (T-test, p=0.735), indicating minor influence on water quality of WWTP operation. However, greater divergences of KALTRIX values were observed after WWTP operation. So, it could be assumed that KALTRIX is more sensitive, and is able to discriminate minor trends not easily detected by TRIX.

The same data set from Kalamitsi and data from South Evoikos used as a reference site (which has been reported as a having high-good ecological quality status, directive 2000/60/EC) were applied to the UNTRIX index (equation 2). Using the box and whisker plots (Table 1), we defined the 100th, 75th and 50th plots of the Southern Evoikos and Kalamitsi areas before and after the WWPT operation.


**Table 1.** The 100th, 75th and 50th box and whisker plots of UNTRIX values from S. Evoikos (reference site) and Kalamitsi before and after the WWTP operation

According to the TQRTRIX index (equation 3), water quality is classified as "moderate" both before (0.56) and after (0.59) the WWTP operation. On the other hand, according to the boxplot procedure [10], water quality is classified as "moderate" before the WWTP operation (75th of the site = 3.142>3.097=100th of the reference), while after the WWTP operation it is classified as "good" [75th of the Kalamitsi =2.994>2.118=75th of the reference (S. Evoikos) and simultaneously 100th of the Kalamitsi = 4.063>3.097=100th of the reference (S. Evoikos)].

#### **2.2. Piecewise linear regression models with breakpoints**

#### *2.2.1. Description of the method*

As the concentrations of Chl-a and DIN from Kalamitsi displayed a profile much lower than those of the Northern Adriatic Sea [11, 7], the TRIX was re-estimated (KALTRIX) according to

**•** zero values of the variables were substituted by the detection limit value of the applied

**•** extreme values (m ± 2.5 SD) were excluded (where m is the mean value and SD is the

KALTRIX values displayed high temporal variability (0.8-9 units) before the WWPT operation, whereas both temporal and spatial variability occurred after the WWPT operation. In both cases, KALTRIX values were higher than those of TRIX, although the same data set was applied. The mean KALTRIX value was 5.3 before WWPT operation and 5.0 after it, with no statistically significant differences (T-test, p=0.735), indicating minor influence on water quality of WWTP operation. However, greater divergences of KALTRIX values were observed after WWTP operation. So, it could be assumed that KALTRIX is more sensitive, and is able

The same data set from Kalamitsi and data from South Evoikos used as a reference site (which has been reported as a having high-good ecological quality status, directive 2000/60/EC) were applied to the UNTRIX index (equation 2). Using the box and whisker plots (Table 1), we defined the 100th, 75th and 50th plots of the Southern Evoikos and Kalamitsi areas before and

100th 3.097 4.091 4.063 75th 2.118 3.142 2.994 50th 1.768 2.673 2.452

**Table 1.** The 100th, 75th and 50th box and whisker plots of UNTRIX values from S. Evoikos (reference site) and Kalamitsi

According to the TQRTRIX index (equation 3), water quality is classified as "moderate" both before (0.56) and after (0.59) the WWTP operation. On the other hand, according to the boxplot procedure [10], water quality is classified as "moderate" before the WWTP operation (75th of the site = 3.142>3.097=100th of the reference), while after the WWTP operation it is classified as "good" [75th of the Kalamitsi =2.994>2.118=75th of the reference (S. Evoikos) and simultaneously 100th of the Kalamitsi = 4.063>3.097=100th of the reference (S. Evoikos)].

**Kalamitsi** (before WWTP operation)

**Kalamitsi** (after WWTP operation)

*KALTRIX = Chl a aD%O DIN TP* 2.5 log \* -\* \* \* - ( ) (0.5) (4)

**•** the range between upper and lower values of each variable was set to 1 log unit

the Kalamitsi data set, following this procedure:

**•** test of normality was performed for log transformations.

to discriminate minor trends not easily detected by TRIX.

**S. Evoikos** (reference site)

standard deviation of the sample)

This procedure furnished equation (4):

after the WWPT operation.

before and after the WWTP operation

analytical method

26 Hydrodynamics - Concepts and Experiments

If an investigation is focused not only on the classification of the coastal waters' quality but also on the assessment of critical thresholds, then piecewise linear regression models with breakpoints may constitute useful tools [19]. Furthermore, piecewise linear regression models with breakpoints are able to express discontinuity appearing in ecosystems, associated either with habitat fragmentations (usually independent variables) or with "decisions" taken by individuals according to their situation [20, 21]. The breakpoint (bt) is a value of either an independent variable (x) or a dependent variable (y) according to the specific characteristics of the ecosystems and the interests of the investigation, where two separate linear regression equations are estimated: one with y values less or equal to bt and the other one with y values greater than bt. The general form is given by equation 5:

$$y = \begin{cases} b\_{01} + b\_{11}\mathbf{x}\_1 + \dots + b\_{n1}\mathbf{x}\_n & \text{for } \mathbf{x} \text{ or } \mathbf{y} \le \mathbf{b}\mathbf{t} \\ b\_{02} + b\_{12}\mathbf{x}\_1 + \dots + b\_{n2}\mathbf{x}\_n & \text{for } \mathbf{x} \text{ or } \mathbf{y} > \mathbf{b}\mathbf{t} \end{cases} \tag{5}$$

Such models have been used, for example, to estimate the production of crops, to search for thresholds in coastal mangrove forests and to examine a limpet feeding rate and its response to temperature [22-24].

#### *2.2.2. Case study*

As chlorophyll-a and nutrients are two of the most common indicators for eutrophication assessment [25], a piecewise linear regression model was applied with data (275 observations) for chlorophyll-a and nitrates (NO3-N) from the Thermaikos Gulf, Greece, which were collected over a period of two years [20].

Nitrogen was chosen because it is used by cells for the biosynthesis of molecules like chloro‐ phyll-a, DNA, RNA and proteins. Thus, it could be assumed that discontinuity might have appeared because each phytoplankton individual "decides" the way that a given amount of nitrogen will be used, depending on the stage of its cell cycle and cell demands. In the fitted model, the concentration of dissolved nitrogen (NO3-N, μg/l) is the independent variable, while the concentration of chlorophyll-a (mg/m3 ) is the dependent one.

$$\text{Chl}-\text{a} = \begin{cases} a\_1 + b\_1 \text{ \* NO}\_3 & \text{if } \text{Chl-a} \le \text{bt} \\ a\_2 + b\_2 \text{ \* NO}\_3 & \text{if } \text{Chl-a} > \text{bt} \end{cases} \tag{6}$$

Because the above model is piecewise linear, the parameters were estimated using the quasi-Newton method, an interactive method that minimizes the least squares loss function [∑ (observed-predicted)2 ] through iterative convergence of the predefined empirical equation.

According to the present data, the breakpoint that minimizes the loss function and simulta‐ neously gives the highest correlation coefficient (r=0.84) is 2.5 mg/m3 of chlorophyll-a. So the model of equation (6) assumes the form:

$$\text{Cbl} - a = \begin{cases} 1.3787 - 0.08441 \text{ \* NO}\_3 & \text{if } \text{Cbl-a} \le 2.5\\ 4.39438 - 0.16438 \text{ \* NO}\_3 & \text{if } \text{Cbl-a} > 2.5 \end{cases}$$

For all the estimated parameters, the t-test was implemented with:

$$\begin{aligned} \mathbf{H}\_0 &\text{parameter} = 0 \\ \mathbf{H}\_1 &\text{parameter} \neq 0 \end{aligned}$$

In all cases, H0 is rejected as p<0.05. The fitted model explains about 71 % of the observed variability (r2 =0.706).

Plotting the predicted and observed values [20], it was found that the predicted values of Chla were overestimated when observed values were lower than 0.9 (mg/m3 ), and underestimated when observed values were higher than 4.4 (mg/m3 ). Thus, we simulated the model starting with a lower breakpoint value (0.5 mg/m3 Chl-a) and terminated it at breakpoint = 4.8 mg/m3 . The results show that when breakpoint values are lower than 0.9 or higher than 4.4 mg/m3 Chla, no relationship could be estimated between Chl-a and NO3. However, when 0.9≤break‐ point≤4.4 mg/m3 Chl-a, a relationship between both parameters appeared. The negative slopes bi (-0.08411 and -0.16438 respectively) indicate a "tendency of the ecosystem to be stabilized" at the values of 2.5 or 0.9 mg/m3 of Chl-a, when its concentration in the Thermaikos Gulf varies between 2.5 and 4.4, and 0.9 and 2.5, respectively. It can be assumed that, when Chl-a values are higher than 2.5, then competition and predation may tend to stabilize the phytoplankton community, while processes like decomposition and inflows may increase NO3 concentration at the same time.

The estimated critical values of breakpoints (0.9, 2.5 and 4.4) are very close to the criteria for trophic assessment in marine and coastal ecosystems given by the Swedish Environmental Protection Agency (<**1.5**, 1.5-2.2, **2.2-3.2**, 3.2-**5** and >5 μg/l Chl-a) and the low risk "trigger values" (2-4 μg/l Chl-a) which were proposed by the Australian and New Zealand Environ‐ ment and Conservation Council [26]. So, it could be assumed that, in the case of the Thermaikos Gulf, values of Chl-a lower than 0.9 may indicate oligotrophic waters, from 0.9 to 2.5, meso‐ trophic, from 2.5 to 4.4, eutrophic, and above 4.4, hypertrophic waters. Furthermore, it has been reported that Chl-a concentrations between 2 and 4.44 appeared during algal blooms [27,28], so it could be suggested that values of Chl-a close to the breakpoint may indicate the start of an algal bloom, although more field and laboratory investigations are needed to confirm this hypothesis.

#### **2.3. A numerical dynamic hydrobiological model**

#### *2.3.1. Description of the method*

In this part of the chapter, a hydrobiological/bio-hydrodynamic model is described. This model is constituted by two parts: a first hydromechanical part, which computes the hydrodynamics and the matter transfer in a coastal geophysical basin; and a second hydrobiological part, which computes the cells' growth rate (generation of new mass) and the cells' decay (loss of mass). When shallow waters characterize a coastal basin, a two-dimensional, depth-averaged, hydrodynamic model is considered to be quite sufficient for the simulation of the seawater circulation. In this case, the hydrobiological model is also a depth-averaged model, like the one presented below.

3 3

1 3787 0 08441 if Chl-a 2.5 4 39438 0 16438 if Chl-a 2.5

*. . \* NO Chl a*

0 1

For all the estimated parameters, the t-test was implemented with:

a were overestimated when observed values were lower than 0.9 (mg/m3

when observed values were higher than 4.4 (mg/m3

variability (r2

point≤4.4 mg/m3

at the same time.

confirm this hypothesis.

*2.3.1. Description of the method*

**2.3. A numerical dynamic hydrobiological model**

bi

=0.706).

28 Hydrodynamics - Concepts and Experiments

at the values of 2.5 or 0.9 mg/m3

*. . \* NO* ìï - £ - = <sup>í</sup> ï - > î

> H : parameter = 0 H : parameter 0 ¹

In all cases, H0 is rejected as p<0.05. The fitted model explains about 71 % of the observed

Plotting the predicted and observed values [20], it was found that the predicted values of Chl-

with a lower breakpoint value (0.5 mg/m3 Chl-a) and terminated it at breakpoint = 4.8 mg/m3

The results show that when breakpoint values are lower than 0.9 or higher than 4.4 mg/m3 Chla, no relationship could be estimated between Chl-a and NO3. However, when 0.9≤break‐

(-0.08411 and -0.16438 respectively) indicate a "tendency of the ecosystem to be stabilized"

between 2.5 and 4.4, and 0.9 and 2.5, respectively. It can be assumed that, when Chl-a values are higher than 2.5, then competition and predation may tend to stabilize the phytoplankton community, while processes like decomposition and inflows may increase NO3 concentration

The estimated critical values of breakpoints (0.9, 2.5 and 4.4) are very close to the criteria for trophic assessment in marine and coastal ecosystems given by the Swedish Environmental Protection Agency (<**1.5**, 1.5-2.2, **2.2-3.2**, 3.2-**5** and >5 μg/l Chl-a) and the low risk "trigger values" (2-4 μg/l Chl-a) which were proposed by the Australian and New Zealand Environ‐ ment and Conservation Council [26]. So, it could be assumed that, in the case of the Thermaikos Gulf, values of Chl-a lower than 0.9 may indicate oligotrophic waters, from 0.9 to 2.5, meso‐ trophic, from 2.5 to 4.4, eutrophic, and above 4.4, hypertrophic waters. Furthermore, it has been reported that Chl-a concentrations between 2 and 4.44 appeared during algal blooms [27,28], so it could be suggested that values of Chl-a close to the breakpoint may indicate the start of an algal bloom, although more field and laboratory investigations are needed to

In this part of the chapter, a hydrobiological/bio-hydrodynamic model is described. This model is constituted by two parts: a first hydromechanical part, which computes the hydrodynamics

Chl-a, a relationship between both parameters appeared. The negative slopes

of Chl-a, when its concentration in the Thermaikos Gulf varies

), and underestimated

.

). Thus, we simulated the model starting

As mentioned above, **the first part of the model refers to the hydromechanical processes.** These hydromechanical processes concern (a) the hydrodynamics, e.g., the computation of the velocity field and the sea surface elevation field; and (b) the matter transport, described by the advection and dispersion processes.

**Hydrodynamics**. The hydrodynamic-hydromechanical model is based on the usual equations for the conservation of mass and momentum, as described in equations (7) [29, 30, 6]:

$$\begin{aligned} \frac{\partial \mathcal{U}}{\partial t} + \mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathbf{x}} + \mathcal{V}\frac{\partial \mathcal{U}}{\partial y} &= -g\frac{\partial \mathcal{L}}{\partial \mathbf{x}} + f\mathcal{V} + \frac{\mathsf{r}\_{sv}}{\rho}\frac{}{h} - \frac{\mathsf{r}\_{bv}}{\rho} + \nu\_h\frac{\partial^2 \mathcal{U}}{\partial \mathbf{x}^2} + \nu\_h\frac{\partial^2 \mathcal{U}}{\partial y^2} &= a\\ \frac{\partial \mathcal{V}}{\partial t} + \mathcal{U}\frac{\partial \mathcal{V}}{\partial \mathbf{x}} + \mathcal{V}\frac{\partial \mathcal{V}}{\partial y} &= -g\frac{\partial \mathcal{L}}{\partial y} - f\mathcal{U} + \frac{\mathsf{r}\_{sy}}{\rho}\frac{}{h} - \frac{\mathsf{r}\_{by}}{\rho} + \nu\_h\frac{\partial^2 \mathcal{V}}{\partial \mathbf{x}^2} + \nu\_h\frac{\partial^2 \mathcal{V}}{\partial y^2} &= b\\ \frac{\partial \mathcal{L}}{\partial t} + \frac{\partial (\mathcal{U}h)}{\partial \mathbf{x}} + \frac{\partial (\mathcal{V}h)}{\partial \mathbf{y}} &= 0 & c \end{aligned} \tag{7}$$

where h is the seawater depth, U and V are the vertically averaged horizontal velocities, ζ is the sea surface elevation, f is the Coriolis parameter, τsx and τsy are the wind surface shear stresses, τbx and τby are the bottom shear stresses, ν<sup>h</sup> is the dispersion coefficient related to the local vorticity, ρ is the density of the water and g is the gravity acceleration [29, 31].

The output from the hydrodynamic simulation, e.g., the velocity field, is then used as an input in the transport simulation. The specific transport model applied here is based on the tracer method which can be described as a Lagrange-Monte Carlo simulation (random walk method) [29, 3]. The simulation of the transport processes based on this method has the important advantage of avoiding numerical problems. Advection and dispersion processes are simulated by means of this method [32-36]. In this specific work, the two parts of the hydromechanical model, i.e., the hydrodynamic and transport processes, are fully coupled, including the biological processes at the same time. In more detail, the transport part of the model uses, as already mentioned, the velocity field, produced by the hydrodynamic model at each time step. The transport model is described analytically in the following paragraphs.

**Matter transport.** According to the tracer method, a large number of particles that represents a particular amount of mass are released in the marine environment through a source. The location where an algal bloom episode occurs determines the position of that source. The transport and fate of each particle is traced with time. Advection of the particulate matter is computed by the local seawater velocity while turbulent diffusion is simulated by the random Brownian motion of the particles. More specifically, the tracer method is described by the following steps [29, 6]:


$$\mathcal{U}\_r = \sqrt{6D\_{\text{eff}}}\tag{8}$$

where D is the local diffusion coefficient and dt is the time step


The horizontal positions of the particles are computed from the superposition of the deter‐ ministic and stochastic displacements:

$$\mathbf{x}\_{i}^{n+1} = \mathbf{x}\_{i}^{n} + \Delta \mathbf{x}\_{i}^{n} + \Delta \mathbf{x}\_{i}^{n'}\\
\text{with } \Delta \mathbf{x}\_{i}^{n} = \mathbf{u}\_{i}^{n} \begin{pmatrix} \mathbf{x}\_{i}^{n} \ \mathbf{t}^{n} \end{pmatrix} \\
\text{dt and} \\
\Delta \mathbf{x}\_{i}^{n} \ \text{"= } \mathbf{u}\_{i}^{n'} \ \text{dt } \operatorname{rnd}[-1, 1] \tag{9}$$

$$\boldsymbol{y}\_{i}^{n+1} = \boldsymbol{y}\_{i}^{n} + \boldsymbol{\Delta}\boldsymbol{y}\_{i}^{n} + \boldsymbol{\Delta}\boldsymbol{y}\_{i}^{n'} \text{ with } \boldsymbol{\Delta}\boldsymbol{y}\_{i}^{n} = \boldsymbol{v}\_{i}^{n}\left(\boldsymbol{y}\_{i}^{n}, \boldsymbol{t}^{n}\right) \text{ } dt \, \boldsymbol{\Delta}\boldsymbol{y}\_{i}^{n} \text{ } \text{' } = \boldsymbol{v}\_{i}^{n'} \text{ } dt \, \operatorname{rnd}[-1, 1] \tag{10}$$

where Δx<sup>i</sup> <sup>n</sup> and Δy<sup>i</sup> <sup>n</sup> are the deterministic displacements and Δx<sup>i</sup> <sup>n</sup>΄ and Δy<sup>i</sup> <sup>n</sup>΄ are the stochastic displacements, with ui <sup>n</sup> (xi n, tn) and vi <sup>n</sup> (yi n, tn) being the deterministic velocities at time tn at the location xi n and yi <sup>n</sup> of the i particle, ui <sup>n</sup>΄& v<sup>i</sup> <sup>n</sup>΄are the random (stochastic) horizontal velocities at time tn at the locations xi and yi , respectively, ui <sup>n</sup>΄= vi <sup>n</sup>΄ = √(6Dh/dt), with Dh being the horizontal particle diffusion coefficient, and rnd is a random variable distributed uniformly between -1 and +1. Eventually, the spatial particle distribution allows for the computation of the particle concentrations in each grid box.

#### **The second part of the model refers to the hydrobiological processes**

If we define


then these three factors regulate the density of population N, which, in a continuous model, follows an exponential growth, according to equation (11):

$$N = N\_0 \cdot e^{(\mu - \text{TL})t} \tag{11}$$

where t is the time (in days), TL is influenced by various processes such as predation, sinking, water movement and stratification, and μ is dependent on abiotic factors such as temperature, light and nutrients, since they act directly on the individuals' biochemistry [37, 38].

**a.** the velocity field is determined by a set of values of the velocity components at specific

**d.** in the case of instantaneous discharge, a specific amount of particles is released at once from the initial source; while in the case of continuous discharge, for each time step a

**e.** integration in time is executed and the new coordinates of each particle are computed; the motion of each particle is analysed considering: i) a deterministic part which concerns the advective transport, and ii) a stochastic part which concerns the transport-spreading due

The horizontal positions of the particles are computed from the superposition of the deter‐

*i i i i i ii i i x x x x x u x t dt and x u dt rnd* <sup>+</sup> = +D +D D = D= - (9)

*i i i i i ii i i y y y y y v y t dt and y v dt rnd* <sup>+</sup> = +D +D D = D= - (10)

<sup>n</sup>΄= vi

( ) 1 ' ' with , ' [ 1,1] *n n n n n nnn n n*

( ) 1 ' ' with , ' [ 1,1] *n n n n n nnn n n*

, respectively, ui

horizontal particle diffusion coefficient, and rnd is a random variable distributed uniformly between -1 and +1. Eventually, the spatial particle distribution allows for the computation of

then these three factors regulate the density of population N, which, in a continuous model,

<sup>n</sup> are the deterministic displacements and Δx<sup>i</sup>

<sup>n</sup> (yi

<sup>n</sup>΄& v<sup>i</sup>

*dt* <sup>=</sup> (8)

<sup>n</sup>΄ and Δy<sup>i</sup>

n, tn) being the deterministic velocities at time tn at the

<sup>n</sup>΄are the random (stochastic) horizontal velocities

<sup>n</sup>΄ = √(6Dh/dt), with Dh being the

<sup>n</sup>΄ are the stochastic

**c.** the range of the random velocity ± Ur is computed from the following equation (8):

<sup>6</sup> *<sup>r</sup> <sup>U</sup> <sup>D</sup>*

grid points computed by the hydrodynamic part of the model

where D is the local diffusion coefficient and dt is the time step

specific number of particles is introduced in the study domain

**b.** a suitable time step is selected

30 Hydrodynamics - Concepts and Experiments

to diffusion processes.

<sup>n</sup> and Δy<sup>i</sup>

n and yi

at time tn at the locations xi

displacements, with ui

**•** TL as the total losses

where Δx<sup>i</sup>

location xi

If we define

ministic and stochastic displacements:

<sup>n</sup> (xi

the particle concentrations in each grid box.

**•** μ as the growth rate (divisions/day) and

<sup>n</sup> of the i particle, ui

n, tn) and vi

**•** N0 as the starting population of a phytoplanktonic species

follows an exponential growth, according to equation (11):

and yi

**The second part of the model refers to the hydrobiological processes**

In the biological part of the model, the maximum growth rate of phytoplankton μ(T) is consid‐ ered to follow Eppley's equation [39, 40]. Then light intensity and sun shine duration are considered to limit μ(T) according to Steele's equation [41]. The sunshine duration and the photosynthetic active radiation (par) were estimated using field data [42, 43]. Nitrogen and phosphorus limitations (Nlim and Plim respectively) are computed according to Michaelis-Menten enzyme Kinetics. Temperature and nutrients concentrations are considered to take random values, following normal distribution which was estimated from field data [27].


The abovementioned processes, with the relevant factors as functions of different elements and the relevant constants, are given in Tables 2a and 2b.

**Table 2.** (a) processes used in the biological compartment of the model (from: [6]) b) constants used in the relations of the previous table (from: [6])

In more detail, the information for the growth rate coefficients was introduced into the model as a time series of mean daily growth rates. That is, it was computed in a discrete way, instead of using the continuous equation (11). The modelling approach, concerning the increase of the particles due to division and the decrease due to particle losses, is described by the following steps [6]:


$$N\_t = N\_0 + \mu \cdot N\_0 - \left[ \text{TL} \cdot \left( N\_0 + \mu \cdot N\_0 \right) \right] \Rightarrow N\_t = N\_0 \left( 1 + \mu - \text{TL} - \mu \cdot \text{TL} \right) \tag{12}$$

Concerning the boundary conditions for the transport model, it should be pointed out that if the particles reach the coast, they return to their previous position, while the particles reaching the open sea boundary are trapped there and no further computation is made for them; these particles are then excluded from the computational loops. The position of the source is determined from the place where an episode of an algal bloom is observed. Particles' concen‐ trations are then computed from the number of particles, counted at each grid box. As far as the initial conditions are concerned, the numerical simulations always start from scratch, i.e., zero current velocities at the starting point of the simulation time.

#### *2.3.2. Case studies*

#### **The hydrobiological model**

The application of the mathematical hydrobiological model described in the previous section is presented here.

**a. Application to a real episode in the Thermaikos Gulf**. The model runs were based on meteorological data for a particular period in which a real episode of algal bloom took place in the Thermaikos Gulf, starting from January 10th, 2000. The Thermaikos Gulf is a semi-enclosed coastal area, located in the north-western corner of the Aegean Sea, in the east of the Mediterranean Sea, as depicted in Figs. 2 and 3. The Axios, Aliakmon and Loudias are the three main rivers that discharge into the western and north-western coasts of the gulf, while the River Pinios outflows much further south, in the south-western coast of the outer gulf. During recent decades, there has been a major reduction in the total freshwater input to the gulf, mainly as the result of the extraction of river water for irrigation, and the construction of a series of hydroelectric power dams on the Aliakmon [55, 56].

In more detail, the information for the growth rate coefficients was introduced into the model as a time series of mean daily growth rates. That is, it was computed in a discrete way, instead of using the continuous equation (11). The modelling approach, concerning the increase of the particles due to division and the decrease due to particle losses, is described by the following

**a.** After a time period of one day, a new number of particles equal to (μ∙N0) are generated. In that product (i.e. μ∙N0) μ is the growth rate coefficient and N0 is the number of particles before the division, denoting the total number of particles of the previous day. The "new" total number of particles N, after the growth process, is afterward computed from the

**b.** The new particles' position is determined from the positions of other particles randomly

**c.** Additionally, a decay coefficient (TL) is adopted for the population decrease and particle disappearance of the water column. So, at the end of a period of a day, a number of particles equal to TL∙N is drawn out, where TL is the aforementioned decay coefficient and N is the total number of particles of the previous day, plus the particles added due to the growth process of the present day, e.g., TL∙N = TL∙(N0+μ∙N0). These particles are no

**d.** Finally, the total number of particles Nt, after growth and decay process, is given by eq. 13:

Concerning the boundary conditions for the transport model, it should be pointed out that if the particles reach the coast, they return to their previous position, while the particles reaching the open sea boundary are trapped there and no further computation is made for them; these particles are then excluded from the computational loops. The position of the source is determined from the place where an episode of an algal bloom is observed. Particles' concen‐ trations are then computed from the number of particles, counted at each grid box. As far as the initial conditions are concerned, the numerical simulations always start from scratch, i.e.,

The application of the mathematical hydrobiological model described in the previous section

**a. Application to a real episode in the Thermaikos Gulf**. The model runs were based on meteorological data for a particular period in which a real episode of algal bloom took place in the Thermaikos Gulf, starting from January 10th, 2000. The Thermaikos Gulf is a semi-enclosed coastal area, located in the north-western corner of the Aegean Sea, in the

m

ë û (12)

 m

00 00 0 (1 ( ) ) *N N N LN N N N L L t t* = + × - T × + × => = + - T - × T

 mé ù

zero current velocities at the starting point of the simulation time.

steps [6]:

equation N=N0 + μ∙N0.

32 Hydrodynamics - Concepts and Experiments

longer taken into account in the computations.

m

selected.

*2.3.2. Case studies*

is presented here.

**The hydrobiological model**

The hydrodynamic model was fully coupled with the transport and biological models described in the previous sections. The equations of the hydrodynamic model were numeri‐ cally solved by the finite difference method on an orthogonal, staggered grid (Arakawa C grid), where U and V components refer to the nodes' sides and ζ refers to the interior of each mesh. The field of the gulf was discretized using 32×36 grid cells. The spatial-discretization step was dx = 1000 m.

**Figure 2.** The Thermaikos Gulf (extending from Thessaloniki's coastline to the upper dash line) and the extended basin of the gulf (extending from Thessaloniki's coastline to the lower dash line) (Google Earth)

Fig. 3a depicts the study area of the Thermaikos basin, with A being the initial location of the appearance of an HAB episode, B the control position-location with reported measured values of HAB concentrations, T the city of Thessaloniki, P the port of Thessaloniki, E the area called Megalo Emvolo, and M1, M2 and M3 the mussel culture areas lying at the west coasts of the gulf. Fig. 3b depicts the grid described above, as the platform for the numerical solution of the equations of the model.

The present application to the real coastal basin of Thermaikos Gulf tried to simulate the fate of harmful algal cells after a sudden appearance of HAB, using *Dinophysis* abundances 104 -105 cells/lit recorded in north-eastern Thermaikos on January 10th, 2000 (position A in Fig. 3a). For the simulation, a number of 1000 particles was initially used at the location (A), and 619,000 particles were finally counted 21 days later (this latter number included particles that may have dispersed in the coastal basin of the gulf or escaped out of the gulf). The model runs led ultimately to the distribution of *Dinophysis* populations in the coastal basin of the gulf. Two

Figure 3. (a) Thermaikos Gulf, where A and B are control positions for the model, T the city of Thessaloniki, P the port of Thessaloniki, E the area called Megalo Emvolo, M1, M2 and M3 mussel culture areas and (b) the grid of the basin's discretization. The axes' coordinates are in number of spatial steps (from: [6]). **Figure 3.** (a) Thermaikos Gulf, where A and B are control positions for the model, T the city of Thessaloniki, P the port of Thessaloniki, E the area called Megalo Emvolo, M1, M2 and M3 mussel culture areas and (b) the grid of the basin's discretization. The axes' coordinates are in number of spatial steps (from: [6]).

control points were used corresponding to the locations A and B, as depicted in Fig. 2. The hydrodynamic circulation of the gulf was activated by variable winds prevailing over this period (between January 10th and 31st, 2000, with a predominance of northern winds that correspond to the most frequent conditions in the Thermaikos Gulf). Concerning the effect of tide on the general pattern of dispersion, it is considered negligible, since the mean tidal signal of 25 cm is quite small and the residual currents due to the tide are negligible. Moreover, riverine water inputs and seawater density differences in general were not taken into account, because the total freshwater input of the rivers has been dramatically reduced in the last decades. Finally, the present numerical simulation was performed for wind-generated circulation, under the winds blowing over the area during the period between January 10th and 31st, 2000. The time series of wind data was recorded from a meteorological station of the Institute of Forestry Research at the location of Sani, Chalkidiki. As mentioned, the transport model was applied considering the position A as the local source of particulate matter and instantaneous appearance of algal bloom, where large concentrations of *Dinophysis* (HAB) were recorded in the field. The present application to the real coastal basin of Thermaikos Gulf tried to simulate the fate of harmful algal cells after a sudden appearance of HAB, using *Dinophysis* abundances 104-105 cells/lit recorded in north-eastern Thermaikos on January 10th, 2000 (position A in Fig. 3a). For the simulation, a number of 1000 particles was initially used at the location (A), and 619,000 particles were finally counted 21 days later (this latter number included particles that may have dispersed in the coastal basin of the gulf or escaped out of the gulf). The model runs led ultimately to the distribution of *Dinophysis* populations in the coastal basin of the gulf. Two control points were used corresponding to the locations A and B, as depicted in Fig. 2. The hydrodynamic circulation of the gulf was activated by variable winds prevailing over this period (between January 10th and 31st, 2000, with a predominance of northern winds that correspond to the most frequent conditions in the Thermaikos Gulf). Concerning the effect of tide on the general pattern of dispersion, it is considered negligible, since the mean tidal signal of 25 cm is quite small and the residual currents due to the tide are negligible. Moreover, riverine water inputs and seawater density differences in general were not taken into account, because the total freshwater input of the rivers has been dramatically reduced in the last decades. Finally, the present numerical simulation was performed for wind-generated circulation, under the winds blowing over the area during the period between January 10th and 31st, 2000. The time series of wind data was recorded from a meteorological station of the Institute of Forestry Research at the location of Sani, Chalkidiki. As mentioned, the transport model was applied

The daily variable estimated growth rate used in the biological model has a range between 0.176 and 0.489 divisions per day [6]. These values are comparable with those estimated using the long-term (one year) deterministic model for dinoflagellates in the same area [37, 38]. In the literature, the values of dinoflagellates' growth rates vary within a range of 0.04-0.7 divisions per day, depending on time and local conditions [1, 48, 57]. considering the position A as the local source of particulate matter and instantaneous appearance of algal bloom, where large concentrations of *Dinophysis* (HAB) were recorded in the field. The daily variable estimated growth rate used in the biological model has a range between 0.176 and 0.489 divisions per day [6]. These values are comparable with those estimated using the long-term (one year) deterministic model for dinoflagellates in the same area [37, 38]. In the literature, the values of dinoflagellates' growth rates vary within a range

The model results concerning the distribution of the harmful phytoplankton cells for time periods of five days, 10 days and 21 days after the most intense bloom of January 10th are given in Figs. 4a, 4b and 4c [6]. According to the abovementioned patterns of the *Dinophysis* dispersion, it may be observed that five days after the bloom at location A (Fig. 4a), *Dinophy‐* of 0.04-0.7 divisions per day, depending on time and local conditions [1, 48, 57]. The model results concerning the distribution of the harmful phytoplankton cells for time periods of five days, 10 days and 21 days after the most intense bloom of January 10th are given in Figs. 4a, 4b and 4c [6]. According to the

are not only aesthetically undesirable but also hazardous on some occasions.

abovementioned patterns of the *Dinophysis* dispersion, it may be observed that five days after the bloom at location A (Fig. 4a), *Dinophysis* masses were transported and concentrated in the area of Megalo Emvolo (position E in Fig. 2) west of position A; 10 days after the initial bloom (Fig. 4b), higher concentrations of *Dinophysis* reached some coastal areas of the west and north-west Thermaikos Gulf, where the largest mussel farms of Greece lie; finally, 21 days after the bloom (Fig. 4c), the present results show that concentrations of *Dinophysis* have been dispersed to a relatively large area of the outer Thermaikos Gulf, and mainly along the west coast of the gulf [6]. Of course, because a large part of the western as well as the eastern coastline consists of bathing coasts, e.g., tourist sea beaches, situations of large algal concentrations *sis* masses were transported and concentrated in the area of Megalo Emvolo (position E in Fig. 2) west of position A; 10 days after the initial bloom (Fig. 4b), higher concentrations of *Dinophysis* reached some coastal areas of the west and north-west Thermaikos Gulf, where the largest mussel farms of Greece lie; finally, 21 days after the bloom (Fig. 4c), the present results show that concentrations of *Dinophysis* have been dispersed to a relatively large area of the outer Thermaikos Gulf, and mainly along the west coast of the gulf [6]. Of course, because a large part of the western as well as the eastern coastline consists of bathing coasts, e.g., tourist sea beaches, situations of large algal concentrations are not only aesthetically undesirable but also hazardous on some occasions.

control points were used corresponding to the locations A and B, as depicted in Fig. 2. The hydrodynamic circulation of the gulf was activated by variable winds prevailing over this period (between January 10th and 31st, 2000, with a predominance of northern winds that correspond to the most frequent conditions in the Thermaikos Gulf). Concerning the effect of tide on the general pattern of dispersion, it is considered negligible, since the mean tidal signal of 25 cm is quite small and the residual currents due to the tide are negligible. Moreover, riverine water inputs and seawater density differences in general were not taken into account, because the total freshwater input of the rivers has been dramatically reduced in the last decades. Finally, the present numerical simulation was performed for wind-generated circulation, under the winds blowing over the area during the period between January 10th and 31st, 2000. The time series of wind data was recorded from a meteorological station of the Institute of Forestry Research at the location of Sani, Chalkidiki. As mentioned, the transport model was applied considering the position A as the local source of particulate matter and instantaneous appearance of algal bloom, where large concentrations of *Dinophysis* (HAB)

**Figure 3.** (a) Thermaikos Gulf, where A and B are control positions for the model, T the city of Thessaloniki, P the port of Thessaloniki, E the area called Megalo Emvolo, M1, M2 and M3 mussel culture areas and (b) the grid of the basin's

5 . 0 0

1 0 . 0 0

1 5 . 0 0

2 0 . 0 0

2 5 . 0 0

3 0 . 0 0 **b**

5 . 0 0 1 0. 0 0 1 5 . 0 0 2 0 . 0 0 2 5 . 0 0 3 0 . 0 0 3 5 . 0 0

Figure 3. (a) Thermaikos Gulf, where A and B are control positions for the model, T the city of Thessaloniki, P the port of Thessaloniki, E the area called Megalo Emvolo, M1, M2 and M3 mussel culture areas and (b) the grid of the basin's discretization. The axes' coordinates

The daily variable estimated growth rate used in the biological model has a range between 0.176 and 0.489 divisions per

abovementioned patterns of the *Dinophysis* dispersion, it may be observed that five days after the bloom at location A (Fig. 4a), *Dinophysis* masses were transported and concentrated in the area of Megalo Emvolo (position E in Fig. 2) west of position A; 10 days after the initial bloom (Fig. 4b), higher concentrations of *Dinophysis* reached some coastal areas of the west and north-west Thermaikos Gulf, where the largest mussel farms of Greece lie; finally, 21 days after the bloom (Fig. 4c), the present results show that concentrations of *Dinophysis* have been dispersed to a relatively large area of the outer Thermaikos Gulf, and mainly along the west coast of the gulf [6]. Of course, because a large part of the western as well as the eastern coastline consists of bathing coasts, e.g., tourist sea beaches, situations of large algal concentrations

**A**

**P T**

The daily variable estimated growth rate used in the biological model has a range between 0.176 and 0.489 divisions per day [6]. These values are comparable with those estimated using the long-term (one year) deterministic model for dinoflagellates in the same area [37, 38]. In the literature, the values of dinoflagellates' growth rates vary within a range of 0.04-0.7

The model results concerning the distribution of the harmful phytoplankton cells for time periods of five days, 10 days and 21 days after the most intense bloom of January 10th are given in Figs. 4a, 4b and 4c [6]. According to the abovementioned patterns of the *Dinophysis* dispersion, it may be observed that five days after the bloom at location A (Fig. 4a), *Dinophy‐*

of 0.04-0.7 divisions per day, depending on time and local conditions [1, 48, 57].

are not only aesthetically undesirable but also hazardous on some occasions.

divisions per day, depending on time and local conditions [1, 48, 57].

large concentrations of *Dinophysis* (HAB) were recorded in the field.

were recorded in the field.

5 10 15 20 25 30 35

discretization. The axes' coordinates are in number of spatial steps (from: [6]).

**E**

5

10

**M3**

Aliacmon

15

20

25

**1 km**

**a**

N

<sup>30</sup> **Thermaikos Gulf**

**M2**

Loudias

**B**

are in number of spatial steps (from: [6]).

Axios

34 Hydrodynamics - Concepts and Experiments

**M1**

Figure 4. The *Dinophysis* concentrations (a) five days, (b) 10 days and (c) 21 days after the bloom in position A (cells/lt). The axes' coordinates are in a number of spatial steps (from: [6]). **Figure 4.** The *Dinophysis* concentrations (a) five days, (b) 10 days and (c) 21 days after the bloom in position A (cells/lt). The axes' coordinates are in a number of spatial steps (from: [6]).

day [6]. These values are comparable with those estimated using the long-term (one year) deterministic model for dinoflagellates in the same area [37, 38]. In the literature, the values of dinoflagellates' growth rates vary within a range The model results concerning the distribution of the harmful phytoplankton cells for time periods of five days, 10 days and 21 days after the most intense bloom of January 10th are given in Figs. 4a, 4b and 4c [6]. According to the In more detail, the simulation shows that 21 days after the bloom at the north-eastern side of the inner Thermaikos Gulf, with an initial concentration of 57,000 cells/lit at position A, and under the influence of prevailing northern and northeastern winds, the main mass was spread to the central and eastern coasts of the inner Thermaikos, as well as to the western coasts of the gulf. Concentrations higher than 200 cells/lt were computed for the region north-east of the Aliacmon river-mouth and south of the Axios river-mouth (position Β in Fig. 2). These concentrations are very close to those reported by [27]. Concerning the north-western coasts of the gulf, e.g., the region of mussel cultures north-east of In more detail, the simulation shows that 21 days after the bloom at the north-eastern side of the inner Thermaikos Gulf, with an initial concentration of 57,000 cells/lit at position A, and under the influence of prevailing northern and north-eastern winds, the main mass was spread to the central and eastern coasts of the inner Thermaikos, as well as to the western coasts of the gulf. Concentrations higher than 200 cells/lt were computed for the region north-east of the Aliacmon river-mouth and south of the Axios river-mouth (position Β in Fig. 2). These

b. Dummy Text **Application of different scenarios to the extended area of the Thermaikos basin**

were not found [6].

step dx = 1852 m (Fig. 5b).

different starting population densities were examined.

the Axios delta (north of position B, area M1), the model resulted to values of up to 400 cells/lt. Furthermore, concerning the west and south-western coasts of the gulf, hosting the mussel cultures south of the Loudias river-mouth (mussel culture area M2, Fig. 2) and west and south of the Aliakmon river-mouth (mussel culture area M3, fig. 2), the of the model resulted to values of up to 900 cells/lt. The *Dinophysis* concentrations reached similar values at the northern area of Megalo Emvolo (position E, fig. 2). Along the coast of Thessaloniki City (position T, Fig. 2), significant *Dinophysis* masses

The application of the aforementioned hydrobiological model, based on different wind conditions prevailing over the extended area of the Thermaikos Gulf, is presented in this section. To this end, the outer part of the external Thermaikos Gulf was also included in the present simulations (Fig. 6a). Furthermore, this application was focused on the investigation of the critical conditions for the appearance of an algal bloom episode and the dispersion of its population in space and time. Additionally, different positions in the Thermaikos Gulf, as sources of the initial population, and

The equations of the hydrodynamic model, in both aforementioned cases, were numerically solved by the finite difference method, on an orthogonal staggered grid (Arakawa C grid) similar to the previous bathymetry (without the greater area of the outer Thermaikos Gulf). The study area was now discretized with a grid of 41×42 cells, with spatial concentrations are very close to those reported by [27]. Concerning the north-western coasts of the gulf, e.g., the region of mussel cultures north-east of the Axios delta (north of position B, area M1), the model resulted to values of up to 400 cells/lt. Furthermore, concerning the west and south-western coasts of the gulf, hosting the mussel cultures south of the Loudias river-mouth (mussel culture area M2, Fig. 2) and west and south of the Aliakmon river-mouth (mussel culture area M3, fig. 2), the of the model resulted to values of up to 900 cells/lt. The *Dinophysis* concentrations reached similar values at the northern area of Megalo Emvolo (position E, fig. 2). Along the coast of Thessaloniki City (position T, Fig. 2), significant *Dinophysis* masses were not found [6].

#### **b. Application of different scenarios to the extended area of the Thermaikos basin**

The application of the aforementioned hydrobiological model, based on different wind conditions prevailing over the extended area of the Thermaikos Gulf, is presented in this section. To this end, the outer part of the external Thermaikos Gulf was also included in the present simulations (Fig. 6a). Furthermore, this application was focused on the investigation of the critical conditions for the appearance of an algal bloom episode and the dispersion of its population in space and time. Additionally, different positions in the Thermaikos Gulf, as sources of the initial population, and different starting population densities were examined.

The equations of the hydrodynamic model, in both aforementioned cases, were numerical‐ ly solved by the finite difference method, on an orthogonal staggered grid (Arakawa C grid) similar to the previous bathymetry (without the greater area of the outer Thermai‐ kos Gulf). The study area was now discretized with a grid of 41×42 cells, with spatial step dx = 1852 m (Fig. 5b).

Figure 5. (a) Thermaikos Gulf, Greece. A and B are the starting points of the simulation, A is the area close to the harbour of Thessaloniki, M1, M2 and M3 are mussel cultures areas [5], and (b) is the grid of the basin's discretization. The axes' coordinates are in number of spatial steps. **Figure 5.** (a) Thermaikos Gulf, Greece. A and B are the starting points of the simulation, A is the area close to the har‐ bour of Thessaloniki, M1, M2 and M3 are mussel cultures areas [5], and (b) is the grid of the basin's discretization. The axes' coordinates are in number of spatial steps.

In this simulation, the following two assumptions drawn from the literature were taken into account:

*Dinophysis spp.* were considered to be the minimum population for starting a bloom episode.

0.000 0.010 0.020 0.030 0.040 0.050 0.060 0.070 0.080 0.090 0.100 0.110 0.120 0.130 0.140 0.150 0.160 0.170 0.180 0.190 0.200

5.00

5

10.00

10

15.00

15

20.00

20

25.00

25

30.00

30

35.00

35

40.00

**b** Hydrodynamic

40

reason, we consider here No = 150 cells/lt.

m/s, are given indicatively in Fig. 6.

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

5 10 15 20 25 30 35 40

5.00

5

number of spatial steps.

10.00

10

15.00

15

20.00

20

25.00

25

30.00

30

35.00

35

40.00

40

**a**

day, while the losses were embodied in the net growth rate.

circulation under north winds 10 m/s

Finally, the combination of the following cases was considered:

a. Although the cell concentration (or abundance) in the peak of the bloom exceeds the number of 104 cells/lt, a minimum of 500 to 1200 cells/lt is a threshold for restrictions in fisheries [44]*.* So, in this study, 2000/lt cells of

b. In almost all the studies, before a bloom episode, *Dinophysis* cell abundance is usually less than 200 cells/lt. For this

The initial model runs' simulation started with a low phytoplankton net growth rate of μ = 0.3 divisions per day. Finally, the simulation tests showed that HAB episodes appeared when the net growth rate takes values close to 1.0 divisions per

A number of 1000 particles was initially used at the location A, and 128,000 particles were finally counted at the end of the simulation. The total simulation time was taken as equal to seven days (time period of one week) which was sufficient for the flow to reach steady state conditions (either for N or for S winds, with speeds of 2 and 10 m/s). The patterns of seawater circulation, corresponding to steady state flow, under the influence of north and south winds of 10

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

Hydrodynamic circulation under south winds 10 m/s

0.000 0.010 0.020 0.030 0.040 0.050 0.060 0.070 0.080 0.090 0.100 0.110 0.120 0.130 0.140 0.150 0.160 0.170 0.180 0.190 0.200 0.210 0.220 0.230 0.240

5 10 15 20 25 30 35 40

Figure 6. Current velocities (m/s) under the influence of (a) N winds (left) and (b) S winds, 10 m/s (right). The axes' coordinates are in

In this simulation, the following two assumptions drawn from the literature were taken into account: Figure 5. (a) Thermaikos Gulf, Greece. A and B are the starting points of the simulation, A is the area close to the harbour of Thessaloniki, M1, M2 and M3 are mussel cultures areas [5], and (b) is the grid of the basin's discretization. The axes' coordinates are in

10.00

15.00

20.00

25.00

30.00

35.00

40.00 **b**

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

**A**

**Thermaikos Gulf**

**1 NM**

**B**

**M1**

**M2 M3**

5.00

10.00

15.00

20.00

25.00

30.00

35.00

40.00

**N**

**a**

concentrations are very close to those reported by [27]. Concerning the north-western coasts of the gulf, e.g., the region of mussel cultures north-east of the Axios delta (north of position B, area M1), the model resulted to values of up to 400 cells/lt. Furthermore, concerning the west and south-western coasts of the gulf, hosting the mussel cultures south of the Loudias river-mouth (mussel culture area M2, Fig. 2) and west and south of the Aliakmon river-mouth (mussel culture area M3, fig. 2), the of the model resulted to values of up to 900 cells/lt. The *Dinophysis* concentrations reached similar values at the northern area of Megalo Emvolo (position E, fig. 2). Along the coast of Thessaloniki City (position T, Fig. 2), significant

**b. Application of different scenarios to the extended area of the Thermaikos basin**

The application of the aforementioned hydrobiological model, based on different wind conditions prevailing over the extended area of the Thermaikos Gulf, is presented in this section. To this end, the outer part of the external Thermaikos Gulf was also included in the present simulations (Fig. 6a). Furthermore, this application was focused on the investigation of the critical conditions for the appearance of an algal bloom episode and the dispersion of its population in space and time. Additionally, different positions in the Thermaikos Gulf, as sources of the initial population, and different starting population densities were examined.

The equations of the hydrodynamic model, in both aforementioned cases, were numerical‐ ly solved by the finite difference method, on an orthogonal staggered grid (Arakawa C grid) similar to the previous bathymetry (without the greater area of the outer Thermai‐ kos Gulf). The study area was now discretized with a grid of 41×42 cells, with spatial step

5.00

**Figure 5.** (a) Thermaikos Gulf, Greece. A and B are the starting points of the simulation, A is the area close to the har‐ bour of Thessaloniki, M1, M2 and M3 are mussel cultures areas [5], and (b) is the grid of the basin's discretization. The

Figure 5. (a) Thermaikos Gulf, Greece. A and B are the starting points of the simulation, A is the area close to the harbour of

a. Although the cell concentration (or abundance) in the peak of the bloom exceeds the number of 104 cells/lt, a minimum of 500 to 1200 cells/lt is a threshold for restrictions in fisheries [44]*.* So, in this study, 2000/lt cells of

b. In almost all the studies, before a bloom episode, *Dinophysis* cell abundance is usually less than 200 cells/lt. For this

The initial model runs' simulation started with a low phytoplankton net growth rate of μ = 0.3 divisions per day. Finally, the simulation tests showed that HAB episodes appeared when the net growth rate takes values close to 1.0 divisions per

A number of 1000 particles was initially used at the location A, and 128,000 particles were finally counted at the end of the simulation. The total simulation time was taken as equal to seven days (time period of one week) which was sufficient for the flow to reach steady state conditions (either for N or for S winds, with speeds of 2 and 10 m/s). The patterns of seawater circulation, corresponding to steady state flow, under the influence of north and south winds of 10

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

Hydrodynamic circulation under south winds 10 m/s

0.000 0.010 0.020 0.030 0.040 0.050 0.060 0.070 0.080 0.090 0.100 0.110 0.120 0.130 0.140 0.150 0.160 0.170 0.180 0.190 0.200 0.210 0.220 0.230 0.240

5 10 15 20 25 30 35 40

Figure 6. Current velocities (m/s) under the influence of (a) N winds (left) and (b) S winds, 10 m/s (right). The axes' coordinates are in

In this simulation, the following two assumptions drawn from the literature were taken into account:

*Dinophysis spp.* were considered to be the minimum population for starting a bloom episode.

0.000 0.010 0.020 0.030 0.040 0.050 0.060 0.070 0.080 0.090 0.100 0.110 0.120 0.130 0.140 0.150 0.160 0.170 0.180 0.190 0.200

5.00

5

10.00

10

15.00

15

20.00

20

25.00

25

30.00

30

35.00

35

40.00

**b** Hydrodynamic

40

10.00

15.00

20.00

25.00

30.00

35.00

40.00 **b**

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

*Dinophysis* masses were not found [6].

36 Hydrodynamics - Concepts and Experiments

dx = 1852 m (Fig. 5b).

**M2 M3**

number of spatial steps.

**A**

**Thermaikos Gulf**

**1 NM**

**B**

**M1**

5.00

5.00

5

number of spatial steps.

10.00

10

15.00

15

20.00

20

25.00

25

30.00

30

35.00

35

40.00

40

**a**

10.00

15.00

20.00

25.00

30.00

35.00

40.00

**N**

**a**

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

axes' coordinates are in number of spatial steps.

reason, we consider here No = 150 cells/lt.

m/s, are given indicatively in Fig. 6.

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

5 10 15 20 25 30 35 40

day, while the losses were embodied in the net growth rate.

circulation under north winds 10 m/s

Finally, the combination of the following cases was considered:


The initial model runs' simulation started with a low phytoplankton net growth rate of μ = 0.3 divisions per day. Finally, the simulation tests showed that HAB episodes appeared when the net growth rate takes values close to 1.0 divisions per day, while the losses were embodied in the net growth rate. b. In almost all the studies, before a bloom episode, *Dinophysis* cell abundance is usually less than 200 cells/lt. For this reason, we consider here No = 150 cells/lt. The initial model runs' simulation started with a low phytoplankton net growth rate of μ = 0.3 divisions per day. Finally, the simulation tests showed that HAB episodes appeared when the net growth rate takes values close to 1.0 divisions per day, while the losses were embodied in the net growth rate.

A number of 1000 particles was initially used at the location A, and 128,000 particles were finally counted at the end of the simulation. The total simulation time was taken as equal to seven days (time period of one week) which was sufficient for the flow to reach steady state conditions (either for N or for S winds, with speeds of 2 and 10 m/s). The patterns of seawater circulation, corresponding to steady state flow, under the influence of north and south winds of 10 m/s, are given indicatively in Fig. 6. A number of 1000 particles was initially used at the location A, and 128,000 particles were finally counted at the end of the simulation. The total simulation time was taken as equal to seven days (time period of one week) which was sufficient for the flow to reach steady state conditions (either for N or for S winds, with speeds of 2 and 10 m/s). The patterns of seawater circulation, corresponding to steady state flow, under the influence of north and south winds of 10 m/s, are given indicatively in Fig. 6.

Figure 6. Current velocities (m/s) under the influence of (a) N winds (left) and (b) S winds, 10 m/s (right). The axes' coordinates are in number of spatial steps. **Figure 6.** Current velocities (m/s) under the influence of (a) N winds (left) and (b) S winds, 10 m/s (right). The axes' coordinates are in number of spatial steps.

Finally, the combination of the following cases was considered: Finally, the combination of the following cases was considered:


study area of the gulf,

**c.** two different wind speeds, the first one corresponding to low winds of 2 m/s and the second one corresponding to relatively strong winds of 10 m/s. c. two different wind speeds, the first one corresponding to low winds of 2 m/s and the second one corresponding to relatively strong winds of 10 m/s. These scenarios led to eight different patterns of algal dispersion, as presented in the following sequence: four patterns

a. two different positions for the starting point, the first one, position A, in the inner part of the gulf, and the second

These scenarios led to eight different patterns of algal dispersion, as presented in the following sequence: four patterns corresponding to the case of initial algal bloom at location A, and four patterns corresponding to the case of initial algal bloom at location B (with A and B depicted in Fig. 5a). corresponding to the case of initial algal bloom at location A, and four patterns corresponding to the case of initial algal bloom at location B (with A and B depicted in Fig. 5a). Dummy Text **Algal bloom at location A**

#### **Algal bloom at location A**

More specifically, in the case of starting an algal bloom at location A in the inner part of the gulf, the patterns of the dispersion of the cells under the influence of north winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 7; the patterns of the dispersion of the cells under the influence of south winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 8. More specifically, in the case of starting an algal bloom at location A in the inner part of the gulf, the patterns of the dispersion of the cells under the influence of north winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 7; the patterns of the dispersion of the cells under the influence of south winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 8.

The following characteristics [5] can be described: The following characteristics [5] can be described:


Figure 7. The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of north winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. **Figure 7.** The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of north winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

#### **Algal bloom at location B**

In the case of starting an algal bloom at location B in the outer part of the gulf, the patterns of the dispersion of the cells under the influence of north winds with wind speeds of 2 and 10 m/

Algae dispersion under south winds 10 m/s

15.00

25.00

30.00

35.00

40.00

Algae dispersion under south winds 2 m/s

 An algal bloom may appear only when wind speed is low (~2 m/s) Figure 8. The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of south winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. **Figure 8.** The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of south winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

s, respectively, are given in Fig. 9; the patterns of the dispersion of the cells under the influence of south winds with wind speeds of 2 and 10 m/s, respectively, are given in Fig. 10. South winds may cause dispersion of phytoplankton cells in the inner gulf. More analytically, when the source is located at position B in the outer gulf, north winds do not seem to cause transport Dummy Text **Algal bloom at location B**

In the case of starting an algal bloom at location B in the outer part of the gulf, the patterns of the dispersion of the cells

The following characteristics [5] can be described: and dispersion of phytoplankton cells at the inner gulf during the period of one week. However, cells can be found along the eastern coasts of the outer gulf. In this case (source at position B), HABs were not observed in the present under the influence of north winds with wind speeds of 2 and 10 m/s, respectively, are given in Fig. 9; the patterns of the

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

5 10 15 20 25 30 35 40

5.00

5

10.00

10

15.00

15

20.00

20

25.00

25

35.00

35

40.00

40

15.00

20.00

25.00

30.00

35.00

40.00

**c.** two different wind speeds, the first one corresponding to low winds of 2 m/s and the

a. two different positions for the starting point, the first one, position A, in the inner part of the gulf, and the second

b. two different wind directions, the north (N) and the south (S) winds, which are generally the most frequent in the

c. two different wind speeds, the first one corresponding to low winds of 2 m/s and the second one corresponding to

A bloom can be started with a low population (150 cells/lt) and a net growth rate close to 1.0 divisions per day, in a

South winds with a low speed of 2 m/s can cause higher cell concentrations than north winds with the same speed.

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

Algae dispersion under north winds 10 m/s

5 10 15 20 25 30 35 40

Figure 7. The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of north winds with speeds of

These scenarios led to eight different patterns of algal dispersion, as presented in the following sequence: four patterns corresponding to the case of initial algal bloom at location A, and four patterns corresponding to the case of initial algal bloom at location B (with A and B depicted

More specifically, in the case of starting an algal bloom at location A in the inner part of the gulf, the patterns of the dispersion of the cells under the influence of north winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 7; the patterns of the dispersion of the cells under the influence of south winds with speeds of 2 and 10 m/s, respectively, are given in Fig. 8.

**•** If the source is located in the Thessaloniki Gulf (inner part of the Thermaikos Gulf), the dispersion of the bloom in the outer part of the Thermaikos Gulf is rather towards the rear **•** A bloom can be started with a low population (150 cells/lt) and a net growth rate close to

the outer part of the Thermaikos Gulf is rather towards the rear

**•** South winds with a low speed of 2 m/s can cause higher cell concentrations than north winds

2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. **Figure 7.** The concentrations of algae (cells/lt) one week after the bloom in position A under the influence of north winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

5.00

In the case of starting an algal bloom at location B in the outer part of the gulf, the patterns of the dispersion of the cells under the influence of north winds with wind speeds of 2 and 10 m/

10.00

15.00

20.00

5101520253035

25.00

30.00

35.00

40.00

40

second one corresponding to relatively strong winds of 10 m/s.

one, position B, in the outer gulf, as depicted in Fig. 5a,

in Fig. 5a).

**Algal bloom at location A**

are given in Fig. 8.

study area of the gulf,

38 Hydrodynamics - Concepts and Experiments

relatively strong winds of 10 m/s.

Dummy Text **Algal bloom at location A**

bloom at location B (with A and B depicted in Fig. 5a).

with the same speed.

5.00

10.00

15.00

20.00

25.00

510152025303540

30.00

35.00

40.00

period of one week

The following characteristics [5] can be described:

The following characteristics [5] can be described:

1.0 divisions per day, in a period of one week

Algae dispersion under north winds 2 m/s

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

5 10 15 20 25 30 35 40

**Algal bloom at location B**


Figure 9. The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of north winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. **Figure 9.** The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of north winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. 30.00 30 2200 2300 2400 2500 2600 30.00 30 220 230 240 250 260

35.00

35

25.00

25

2700 2800

2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

5.00

5

10.00

10

15.00

15

20.00

20

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

270 280

5 10 15 20 25 30 35 40

Figure 9. The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of north winds with speeds of

More analytically, when the source is located at position B in the outer gulf, north winds do not seem to cause transport and dispersion of phytoplankton cells at the inner gulf during the period of one week. However, cells can be found along the eastern coasts of the outer gulf. In this case (source at position B), HABs were not observed in the present calculations, except under the influence of low winds (~2 m/s). South winds may transfer and disperse phyto‐ plankton cells in the inner gulf and close to the areas of mussel culture after a period of one week, but phytoplankton concentrations in all cases are shown to be lower than 1000 cells/lt.

Figure 10. The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of south winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps. **Figure 10.** The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of south winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

#### **2.4. A stochastic markovian model 2.4. A stochastic markovian model**

From

. . .

1 n1(t)

2 n2(t)

k nk(t)

. . .

#### **2.4.1. Model description** *2.4.1. Model description*

As randomness and uncertainty appear in natural ecosystems, Markovian models, as a part of stochastic processes, have been applied in many ecological processes, including succession, population dynamics, energy flow and diversity [37, 58-61]. In Markovian models, the definition of "states" is a critical point and depends on the specific characteristics of the system to be modelled. A state consists of individuals, variables or number of species, while all the possible states of the system correspond to the "state space" C (C={1,2,3...,k}). The vector n(t)=[n1(t), n2(t), ..., nk(t)] gives the expected number of the members on each one of the k "states", for each t (Table 1). The so-called transition probability (pij(t)) is the probability of a member moving from state *i*, at time *t*, to state *j* at time *t+1*. Transition probabilities can be estimated as: Bez Broja ( 1) where, i, j 1,2,3,...k ( ) *ij ij i n t p n t* As randomness and uncertainty appear in natural ecosystems, Markovian models, as a part of stochastic processes, have been applied in many ecological processes, including succession, population dynamics, energy flow and diversity [37, 58-61]. In Markovian models, the definition of "states" is a critical point and depends on the specific characteristics of the system to be modelled. A state consists of individuals, variables or number of species, while all the possible states of the system correspond to the "state space" C (C={1,2,3...,k}). The vector n(t)=[n1(t), n2(t),..., nk(t)] gives the expected number of the members on each one of the k "states", for each t (Table 1). The so-called transition probability (pij(t)) is the probability of a member moving from state *i*, at time *t*, to state *j* at time *t+1*. Transition probabilities can be estimated as:

$$p\_{ij} = \frac{n\_{ij}(t+1)}{n\_i(t)} \\ \text{where, i, j = 1, 2, 3, \dots k}$$

n11(t+1) 11

n21(t+1) 21

> . . .

nk1(t+1) 1

1 ( 1) ( ) *n t*

2 ( 1) ( ) *n t*

*n t*

( 1) ( ) *<sup>k</sup>*

*k n t*

*n t*

*n t*

n12(t+1) 12

n22(t+1) 22

> . . .

nk2(t+1) 2

1 ( 1) ( ) *n t*

2 ( 1) ( ) *n t*

*n t*

( 1) ( ) *<sup>k</sup>*

*k n t*

*n t*

*n t*

. .

. .

. . . . . .

. . .

n1k(t+1) 1

n2k(t+1) 2

> . . .

nkk(t+1)

*n t*

*n t*

1 ( 1) ( ) *<sup>k</sup>*

2 ( 1) ( ) *<sup>k</sup>*

*n t*

 <sup>1</sup> *kk*

*k n t*

*n t*

*n t*

1

2

*kk*

*p*

*p*

*k*

*p*

*k*

12

22

2

*k*

*p*

The "state space" and time can be either continuous or discrete [62, 63]. A Markov process is called stationary or homogeneous if P(t)=pij, *t* , otherwise it is called non-homogeneous. A non-homogeneous Markov process undergoes cyclic behaviour if there exists *d* such that P(ad+s)=P(s), *t* and s={1,2,...,d-1} [64]. For ecosystems with seasonal cycles d=4, let P(0), P(1), P(2) and P(3) be the transition matrices from winter to spring, spring to summer,

*p*

*p*

11

21

1

*k*

*p*

Table 3. Transition probabilities matrix (P) with transition probabilities pij (from: [67])

*p*

*p*


**Table 3.** Transition probabilities matrix (P) with transition probabilities pij (from: [67])

The "state space" and time can be either continuous or discrete [62, 63]. A Markov process is called stationary or homogeneous if P(t)=pij, ∀ *t* ∈ℕ, otherwise it is called non-homogeneous. A non-homogeneous Markov process undergoes cyclic behaviour if there exists *d* ∈ℕ such that P(ad+s)=P(s), ∀ *t* ∈ℕ and s={1,2,...,d-1} [64]. For ecosystems with seasonal cycles d=4, let P(0), P(1), P(2) and P(3) be the transition matrices from winter to spring, spring to summer, summer to autumn and autumn to winter, respectively. Then, the transition matrices from winter to winter (P0), spring to spring (P1), summer to summer (P2) and autumn to autumn (P3) are given by:

> ( ) () ( ) ( ) () ( ) ( ) ( ) ( ) ( ) ( ) () ( ) ( ) () ( ) 0 1 2 3 P = P 0 P 1 P 2 P 3 P = P 1 P 2 P 3 P 0 P = P 2 P 3 P 0 P 1 P = P 3 P 0 P 1 P 2

The matrices Pi , i={0,1,2,3} are stochastic and converge if they are regular, and the convergence is geometrically fast. Their limits (*Pi t* ) are independent of the starting season, meaning that lim(winter to winter)=lim(spring to winter)=... =lim(autumn to winter), and give the probabil‐ ities of the system in "statistical equilibrium" for each season [65-67].

#### *2.4.2. Case study*

More analytically, when the source is located at position B in the outer gulf, north winds do not seem to cause transport and dispersion of phytoplankton cells at the inner gulf during the period of one week. However, cells can be found along the eastern coasts of the outer gulf. In this case (source at position B), HABs were not observed in the present calculations, except under the influence of low winds (~2 m/s). South winds may transfer and disperse phyto‐ plankton cells in the inner gulf and close to the areas of mussel culture after a period of one week, but phytoplankton concentrations in all cases are shown to be lower than 1000 cells/lt.

Figure 10. The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of south winds with speeds of

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

5 10 15 20 25 30 35 40

51015202530354010

Algae dispersion under south winds 10 m/s

As randomness and uncertainty appear in natural ecosystems, Markovian models, as a part of stochastic processes, have been applied in many ecological processes, including succession, population dynamics, energy flow and diversity [37, 58-61]. In Markovian models, the definition of "states" is a critical point and depends on the specific characteristics of the system to be modelled. A state consists of individuals, variables or number of species, while all the possible states of the system correspond to the "state space" C (C={1,2,3...,k}). The vector n(t)=[n1(t), n2(t), ..., nk(t)] gives the expected number of the members on each one of the k "states", for each t (Table 1). The so-called transition probability (pij(t)) is the probability of a member moving from state *i*, at time *t*, to state *j* at time *t+1*. Transition probabilities can be estimated as:

**To states To states at time (t+1)**

( 1) where, i, j 1,2,3,...k ( )

n11(t+1) 11

n21(t+1) 21

> . . .

nk1(t+1) 1

1 ( 1) ( ) *n t*

2 ( 1) ( ) *n t*

*n t*

( 1) ( ) *<sup>k</sup>*

*k n t*

*n t*

*n t*

11

<sup>+</sup> = = *ij*

21

1

*k*

*p*

Table 3. Transition probabilities matrix (P) with transition probabilities pij (from: [67])

*p*

*p*

*i*

*n t*

*n t*

states time (t) 1 2 . . k

n12(t+1) 12

n22(t+1) 22

> . . .

nk2(t+1) 2

1 ( 1) ( ) *n t*

2 ( 1) ( ) *n t*

*n t*

( 1) ( ) *<sup>k</sup>*

*k n t*

*n t*

*n t*

. .

. .

. . . . . .

. . .

n1k(t+1) 1

n2k(t+1) 2

> . . .

nkk(t+1)

*n t*

*n t*

1 ( 1) ( ) *<sup>k</sup>*

2 ( 1) ( ) *<sup>k</sup>*

*n t*

 <sup>1</sup> *kk*

*k n t*

*n t*

*n t*

1

2

*kk*

*p*

*p*

*k*

*p*

*k*

12

22

2

*k*

*p*

The "state space" and time can be either continuous or discrete [62, 63]. A Markov process is called stationary or homogeneous if P(t)=pij, *t* , otherwise it is called non-homogeneous. A non-homogeneous Markov process undergoes cyclic behaviour if there exists *d* such that P(ad+s)=P(s), *t* and s={1,2,...,d-1} [64]. For ecosystems with seasonal cycles d=4, let P(0), P(1), P(2) and P(3) be the transition matrices from winter to spring, spring to summer,

*p*

*p*

2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

5.00

**Figure 10.** The concentrations of algae (cells/lt) one week after the bloom in position B under the influence of south winds with speeds of 2 m/s (left) and 10 m/s (right) (from: [5]). The axes' coordinates are in number of spatial steps.

As randomness and uncertainty appear in natural ecosystems, Markovian models, as a part of stochastic processes, have been applied in many ecological processes, including succession, population dynamics, energy flow and diversity [37, 58-61]. In Markovian models, the definition of "states" is a critical point and depends on the specific characteristics of the system to be modelled. A state consists of individuals, variables or number of species, while all the possible states of the system correspond to the "state space" C (C={1,2,3...,k}). The vector n(t)=[n1(t), n2(t),..., nk(t)] gives the expected number of the members on each one of the k "states", for each t (Table 1). The so-called transition probability (pij(t)) is the probability of a member moving from state *i*, at time *t*, to state *j* at time *t+1*. Transition probabilities can be

10.00

15.00

20.00

25.00

30.00

35.00

40.00

**2.4. A stochastic markovian model**

5.00 10.00 15.00 20.00 25.00 30.00 35.00 40.00

25303540200

40 Hydrodynamics - Concepts and Experiments

Algae dispersion under south winds 2 m/s

5 10 15 20 25 30 35 40

Bez Broja ( 1) where, i, j 1,2,3,...k ( ) *ij*

From

. . .

1 n1(t)

*ij*

*p*

2 n2(t)

k nk(t)

. . .

*i n t*

*n t*

**2.4.1. Model description**

*2.4.1. Model description*

**2.4. A stochastic markovian model**

5.00

5

10.00

15.00

101520

20.00

25.00

30.00

35.00

40.00

*ij*

*p*

estimated as:

A study of zoobenthos dynamics from the Thermaikos Gulf, Greece, is presented in this section [37, 67]. In this study, the following aspects are examined: the seasonal variation of zoobenthos richness (number of species/0.2 m2 ), abundance (individuals/0.2 m2 ) and diversity. Concerning zoobenthos richness, four states were chosen: {A1(up to 10 species), A2 (from 11 to 15 species), A3 (from 16 to 20 species) and A4 (more than 20 species). The definition of the "state space" was based on the literature. For example, fewer than 10 species have been reported in pollutant waters with oxygen concentrations less than 1.4 ml/l [68, 69], while in gulfs similar to the Thermaikos Gulf, without pollution effects, more than 20 species have been reported [70]. The definition of A3 and A4 was based on data from areas under eutrophication processes, where the number of species varied seasonally between 10 and 20 [70-72]. On abundance, four states were defined {(B1 (less than 90 individuals), B2 (from 91-140), B3 (from 141 to 200) and B4 (more than 200)}, following the same assumptions as before when considering diversity. The estimated matrices with transition probabilities from each season to the next, starting from winter to spring, are presented in Table 4.


**Table 4.** Transition probabilities for species richness (a) and abundance (b)

The matrices Pi , i={0, 1, 2, 3} and their limits (*Pi t* ) were estimated (Fig. 11). Additionally, in order to test if the number of individuals (abundance) and the number of species (richness) are independent, correlation coefficients were estimated and, furthermore, simple regression models were applied to each set of data separately. In all cases, the statistical tests confirm that the two variables are independent (p<0.05). So, the joint probability P(Ai and Bi ) in statistical equilibrium could be estimated as follows:

zoobenthos richness, four states were chosen: {A1(up to 10 species), A2 (from 11 to 15 species), A3 (from 16 to 20 species) and A4 (more than 20 species). The definition of the "state space" was based on the literature. For example, fewer than 10 species have been reported in pollutant waters with oxygen concentrations less than 1.4 ml/l [68, 69], while in gulfs similar to the Thermaikos Gulf, without pollution effects, more than 20 species have been reported [70]. The definition of A3 and A4 was based on data from areas under eutrophication processes, where the number of species varied seasonally between 10 and 20 [70-72]. On abundance, four states were defined {(B1 (less than 90 individuals), B2 (from 91-140), B3 (from 141 to 200) and B4 (more than 200)}, following the same assumptions as before when considering diversity. The estimated matrices with transition probabilities from each season to the next, starting from

> **From winter to spring From spring to summer** A1 A2 A3 A4 A1 A2 A3 A4

**From summer to autumn From autumn to winter** A1 A2 A3 A4 A1 A2 A3 A4

**From winter to spring From spring to summer** B1 B2 B3 B4 B1 B2 B3 B4

**From summer to autumn From autumn to winter** B1 B2 B3 B4 B1 B2 B3 B4

A1 0.5 0.5 0 0 0.584 0.25 0.083 0.083 A2 0.714 0.143 0.143 0 0.112 0.444 0.333 0.111 A3 0.461 0.461 0.078 0 0 0.333 0.667 0 A4 0 0.5 0.25 0.25 0 0 0.5 0.5

A1 0.5 0.5 0 0 0.584 0.25 0.083 0.083 A2 0.714 0.143 0.143 0 0.112 0.444 0.333 0.111 A3 0.461 0.461 0.078 0 0 0.333 0.667 0 A4 0 0.5 0.25 0.25 0 0 0.5 0.5 (a)

B1 0.833 0.167 0 0 0.786 0.071 0.143 0 B2 0.712 0.224 0.064 0 0.4 0.2 0.2 0.2 B3 0.222 0.222 0.334 0.222 0.2 0.4 0.2 0.2 B4 0.571 0.143 0.286 0 0 0 0.5 0.5

B1 0.571 0.143 0.214 0.072 0.333 0.167 0.25 0.25 B2 0.4 0.2 0.2 0.2 0.25 0.25 0.25 0.25 B3 0.25 0.25 0.25 0.25 0 0.4 0.4 0.2 B4 0 0 0.333 0.667 0 0 0.6 0.4 (b)

**Table 4.** Transition probabilities for species richness (a) and abundance (b)

winter to spring, are presented in Table 4.

42 Hydrodynamics - Concepts and Experiments

$$\mathbf{P(A\_i \text{and } B\_i) = P(A\_i) \text{ \* } P(B\_i) \to i = \{1, 2, 3, 4\}}.$$

P(Ai and Bi) = P(Ai) \* P(Bi) i={1,2,3,4}.

estimated as follows:

Figure 11. The limits of transition probabilities' matrices for each season: a) for species richness and b) for abundance (population density). The values of the probabilities are represented in the *y* axes. **Figure 11.** The limits of transition probabilities' matrices for each season: a) for species richness and b) for abundance (population density). The values of the probabilities are represented in the *y* axes.

According to Fig. 11, it could be said that winter is the best season, because the probabilities of the system being in the more "rich states" (A3, A4 and B3 and B4) are 0.6916 and 0.6221, respectively. At the same time, the probability of the system to have more than 16 species and more than 141 individuals {(A3 and B3) or (A3 and B4) or (A4 and B3) or (A4 and B4)}, in steady state, is higher than 0.42 (Fig. 12). On the other hand, during spring, the probabilities of A1 (=0.46) and B1 (=0.5221) According to Fig. 11, it could be said that winter is the best season, because the probabilities of the system being in the more "rich states" (A3, A4 and B3 and B4) are 0.6916 and 0.6221, respectively. At the same time, the probability of the system to have more than 16 species and

18

could be estimated as follows:

The matrices Pi, i={0, 1, 2, 3} and their limits (ܲ

Bez Broja P A and B = P A \* P B i= 1,2,3,4 . ii i i

more than 141 individuals {(A3 and B3) or (A3 and B4) or (A4 and B3) or (A4 and B4)}, in steady state, is higher than 0.42 (Fig. 12). On the other hand, during spring, the probabilities of A1 (=0.46) and B1 (=0.5221) taking their highest values among all the seasons (Fig. 11) and, furthermore, the probability of the system having fewer than 15 species and simultaneously less than 90 individuals, is higher than 0.6 (Fig. 12). In relation to these values, it could be suggested that in spring, strong environmental influences tend to reduce both the species richness and the populations density. During summer, it can be observed (Table 3) that transition probabilities to A1 and A2 still have high values, but by contrast transition proba‐ bilities to A3 and A4 are higher than the respective ones in spring. Concerning abundance, it could be observed (Fig. 11) that B1 is dominant, with the highest probability value (=0.5285). The probability of the system having fewer than 15 species and fewer than 90 individuals is 0.44. So, during summer, population density seems to stabilize in low densities, similar to the behaviour observed in the spring. During autumn, the "richest state", A4, takes its highest probability value (=0.1823) and simultaneously A3 and A4 are higher than those of spring and summer. In addition, B1 probability is lower than that observed for summer, although it continues to dominate (Fig. 11). The joint probability of the system having 11-15 species and up to 140 individuals is estimated at 0.353 (Fig. 12). So, it could be considered that a tendency towards recovery starts from autumn and continues to winter. Figure 11. The limits of transition probabilities' matrices for each season: a) for species richness and b) for abundance (population density). The values of the probabilities are represented in the *y* axes. According to Fig. 11, it could be said that winter is the best season, because the probabilities of the system being in the more "rich states" (A3, A4 and B3 and B4) are 0.6916 and 0.6221, respectively. At the same time, the probability of the system to have more than 16 species and more than 141 individuals {(A3 and B3) or (A3 and B4) or (A4 and B3) or (A4 and B4)}, in steady state, is higher than 0.42 (Fig. 12). On the other hand, during spring, the probabilities of A1 (=0.46) and B1 (=0.5221) taking their highest values among all the seasons (Fig. 11) and, furthermore, the probability of the system having fewer than 15 species and simultaneously less than 90 individuals, is higher than 0.6 (Fig. 12). In relation to these values, it could be suggested that in spring, strong environmental influences tend to reduce both the species richness and the populations density. During summer, it can be observed (Table 3) that transition probabilities to A1 and A2 still have high values, but by contrast transition probabilities to A3 and A4 are higher than the respective ones in spring. Concerning abundance, it could be observed (Fig. 11) that B1 is dominant, with the highest probability value (=0.5285). The probability of the system having fewer than 15 species and fewer than 90 individuals is 0.44. So, during summer, population density seems to stabilize in low densities, similar to the behaviour observed in the spring. During autumn, the "richest state", A4, takes its highest probability value (=0.1823) and simultaneously A3 and A4 are higher than those of spring and summer. In addition, B1 probability is lower than that observed for summer, although it continues to dominate (Fig. 11). The joint probability of the system having 11-15 species and up to 140 individuals is estimated at 0.353 (Fig. 12). So, it could be considered that a tendency towards recovery starts from autumn and continues to winter.

௧

individuals (abundance) and the number of species (richness) are independent, correlation coefficients were estimated and, furthermore, simple regression models were applied to each set of data separately. In all cases, the statistical tests confirm that the two variables are independent (p<0.05). So, the joint probability P(Ai and Bi) in statistical equilibrium

) were estimated (Fig. 11). Additionally, in order to test if the number of

**Figure 12.** Joint probabilities Ai and Bi i={1,2,3,4}, for each season separately. The values of the probabilities are repre‐ sented in the *y* axes.

According to Fig. 11, it could be said that winter is the best season, because the probabilities of the system being in the more "rich states" (A3, A4 and B3 and B4) are 0.6916 and 0.6221, respectively. At the same time, the probability of the system to have more than 16 species and more than 141 individuals {(A3 and B3) or (A3 and B4) or (A4 and B3) or (A4 and (=0.5221) taking their highest values among all the seasons (Fig. 11) and, furthermore, the probability of the system having fewer than 15 species and simultaneously less than 90 individuals, is higher than 0.6 (Fig. 12). In relation to these values, it could be suggested that in spring, strong environmental influences tend to reduce both the species richness and the populations density. During summer, it can be observed (Table 3) that transition probabilities to A1 and A2 still have The above discussion points out that the composition of the benthic community of zoobenthos from the Thermaikos Gulf may be regulated by seasonal sequences of extinction and recolo‐ nization. The two "forcing variables" in the Thermaikos Gulf are oxygen concentration and temperature. So, species with tolerance to high temperature, low oxygen and high pollutant concentrations may regulate their births from spring to summer, while species with no tolerance to low oxygen or/and high temperatures may regulate their births from autumn to winter, when the waters are saturated with oxygen. The low values of abundance (population density) were observed during spring to autumn, and may be caused by mass mortalities due to oxygen depletion in the bottom waters [67].

#### The probability of the system having fewer than 15 species and fewer than 90 individuals is 0.44. So, during summer, **3. Conclusions**

more than 141 individuals {(A3 and B3) or (A3 and B4) or (A4 and B3) or (A4 and B4)}, in steady state, is higher than 0.42 (Fig. 12). On the other hand, during spring, the probabilities of A1 (=0.46) and B1 (=0.5221) taking their highest values among all the seasons (Fig. 11) and, furthermore, the probability of the system having fewer than 15 species and simultaneously less than 90 individuals, is higher than 0.6 (Fig. 12). In relation to these values, it could be suggested that in spring, strong environmental influences tend to reduce both the species richness and the populations density. During summer, it can be observed (Table 3) that transition probabilities to A1 and A2 still have high values, but by contrast transition proba‐ bilities to A3 and A4 are higher than the respective ones in spring. Concerning abundance, it could be observed (Fig. 11) that B1 is dominant, with the highest probability value (=0.5285). The probability of the system having fewer than 15 species and fewer than 90 individuals is 0.44. So, during summer, population density seems to stabilize in low densities, similar to the behaviour observed in the spring. During autumn, the "richest state", A4, takes its highest probability value (=0.1823) and simultaneously A3 and A4 are higher than those of spring and summer. In addition, B1 probability is lower than that observed for summer, although it continues to dominate (Fig. 11). The joint probability of the system having 11-15 species and up to 140 individuals is estimated at 0.353 (Fig. 12). So, it could be considered that a tendency

௧

individuals (abundance) and the number of species (richness) are independent, correlation coefficients were estimated and, furthermore, simple regression models were applied to each set of data separately. In all cases, the statistical tests confirm that the two variables are independent (p<0.05). So, the joint probability P(Ai and Bi) in statistical equilibrium

Figure 11. The limits of transition probabilities' matrices for each season: a) for species richness and b) for abundance (population

B4)}, in steady state, is higher than 0.42 (Fig. 12). On the other hand, during spring, the probabilities of A1 (=0.46) and B1

high values, but by contrast transition probabilities to A3 and A4 are higher than the respective ones in spring. Concerning abundance, it could be observed (Fig. 11) that B1 is dominant, with the highest probability value (=0.5285).

population density seems to stabilize in low densities, similar to the behaviour observed in the spring. During autumn,

0.353 (Fig. 12). So, it could be considered that a tendency towards recovery starts from autumn and continues to winter.

A1 A2 A3 A4

A1 A2 A3 A4

i={1,2,3,4}, for each season separately. The values of the probabilities are repre‐

B1 B2 B3 B4

0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 0.22 0.24

> 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2

) were estimated (Fig. 11). Additionally, in order to test if the number of

towards recovery starts from autumn and continues to winter.

Winter Spring

A1 A2 A3 A4

A1 A2 A3 A4

and Bi

Summer Autumn

0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2

> 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2

sented in the *y* axes.

**Figure 12.** Joint probabilities Ai

The matrices Pi, i={0, 1, 2, 3} and their limits (ܲ

Bez Broja P A and B = P A \* P B i= 1,2,3,4 . ii i i

density). The values of the probabilities are represented in the *y* axes.

could be estimated as follows:

44 Hydrodynamics - Concepts and Experiments

the "richest state", A4, takes its highest probability value (=0.1823) and simultaneously A3 and A4 are higher than those of spring and summer. In addition, B1 probability is lower than that observed for summer, although it continues to dominate (Fig. 11). The joint probability of the system having 11-15 species and up to 140 individuals is estimated at The fact that the TRIX index has been applied in many areas suggests its use as a tool not only for the classification of marine waters but also for the comparison between data from different regions. The problem of its calibration arises when marine ecological processes of the study region differ from those of the Northern Adriatic Sea. In such cases, a re-scaled index like the KALTRIX seems to be more sensitive than TRIX to minor changes in the parameters, but as the index was re-scaled its results are not directly comparable with those of other regions. UNTRIX is also a sensitive index, but the problem of the reference site and its characteristics is evident. So, it could be said that if the interest of a manager is to compare a marine ecosystem before and after a strong human influence (like WWTP), then more sensitive indices like KALTRIX or UNTRIX should be used; otherwise, the use of the TRIX index gives the advantage of enabling direct comparisons between different marine ecosystems.

> An alternative approach for marine water classification and management is the use of piecewise linear regression models with breakpoints. Such models are able to express discon‐ tinuities in ecosystems and to detect critical thresholds in coastal ecosystems.

> Concerning the case of the numerical dynamic hydrobiological model, the dispersion of harmful cells after a HAB episode in a coastal basin was investigated. As expected, the wind speed and direction significantly affect the dispersion of an algal bloom. According to the results of the present work, the location of an initial outbreak in combination with the topography of a coastal basin play an important role in the appearance and the dispersion of an algal bloom, because winds with the same speed induce two different patterns of dispersion and different phytoplankton abundances (or concentrations) in space and time. The applica‐ tion of the dynamic hydrobiological model simulating a real episode showed that after a phytoplankton bloom at the north-east coasts of the Thermaikos Gulf (with variable wind forcing over the coastal basin and prevailing northern winds), concentrations of harmful algal cells reached (a) the north-west and west coasts of the gulf where the largest part of the mussel culture of Greece lies, and (b) the west coasts of the gulf with important bathing beaches, within a period of 10 days. For 20 days after the bloom, the mass of *Dinophysis* (HAB) dispersed in quite large areas of the outer Thermaikos Gulf. Concerning the application of the above model to the larger area of the extended Thermaikos Gulf, the model runs showed that, after an algal bloom, the stronger the winds are, the lower the concentration of the dispersed phytoplank‐

tonic cells is expected to be. This means that stronger winds more effectively disperse the masses, "diluting" the large hazardous populations. Moreover, when the source is located in the inner part of the gulf, the dispersion of the bloom does not seem to reach to the outer parts of the basin under north or south winds. From the applications of the dynamic hydrobiological models described here, it can be said that the spreading of phytoplankton cells in the large area of a gulf or a coastal basin can be investigated analytically, so that the essential measures of prevention can be applied for each different case.

Finally, if time-dependent uncertainty and randomness affect the processes under investiga‐ tion, and seasonal cycles are displayed, then non-homogeneous Markovian models with cyclic behaviour may be useful tools for their description. Furthermore, these type of models may be used as an alternative approach to the study of diversity and stability of such marine ecosystems.

## **Acknowledgements**

We owe many thanks to the Forestry Research Institute for providing us with wind data for the period of January 2000. We also wish to thank Dr O. Gotsis-Skretas, Dr Ch. Kontoyiannis and Dr A Pavlidou from the Hellenic Center for Marine Research for providing the data from S. Evoikos Gulf.

### **Author details**

Dimitrios P. Patoucheas1\* and Yiannis G. Savvidis2

\*Address all correspondence to: dpatouh@bio.auth.gr

1 Ministry of Education and Religious Affairs, Secondary Education of Thessaloniki, Thessaloniki, Greece

2 Department of Civil Engineering of Technological Education, Alexander Technological Educational Institute of Thessaloniki, Thessaloniki, Greece

#### **References**

[1] Xie H., Lazure P., Gentien P. Small scale retentive structures and Dinophysis. Journal of Marine Systems 2007: 64(1-4), 173-188.

[2] Ganoulis J. G. Water quality assessment and protection measures of a semi-enclosed coastal area: the Bay of Thermaikos (NE Mediterranean Sea). Mar. Poll. Bull. 1991: 23, 83-87.

tonic cells is expected to be. This means that stronger winds more effectively disperse the masses, "diluting" the large hazardous populations. Moreover, when the source is located in the inner part of the gulf, the dispersion of the bloom does not seem to reach to the outer parts of the basin under north or south winds. From the applications of the dynamic hydrobiological models described here, it can be said that the spreading of phytoplankton cells in the large area of a gulf or a coastal basin can be investigated analytically, so that the essential measures

Finally, if time-dependent uncertainty and randomness affect the processes under investiga‐ tion, and seasonal cycles are displayed, then non-homogeneous Markovian models with cyclic behaviour may be useful tools for their description. Furthermore, these type of models may be used as an alternative approach to the study of diversity and stability of such marine

We owe many thanks to the Forestry Research Institute for providing us with wind data for the period of January 2000. We also wish to thank Dr O. Gotsis-Skretas, Dr Ch. Kontoyiannis and Dr A Pavlidou from the Hellenic Center for Marine Research for providing the data from

1 Ministry of Education and Religious Affairs, Secondary Education of Thessaloniki,

2 Department of Civil Engineering of Technological Education, Alexander Technological

[1] Xie H., Lazure P., Gentien P. Small scale retentive structures and Dinophysis. Journal

of prevention can be applied for each different case.

Dimitrios P. Patoucheas1\* and Yiannis G. Savvidis2

\*Address all correspondence to: dpatouh@bio.auth.gr

Educational Institute of Thessaloniki, Thessaloniki, Greece

of Marine Systems 2007: 64(1-4), 173-188.

ecosystems.

S. Evoikos Gulf.

**Author details**

Thessaloniki, Greece

**References**

**Acknowledgements**

46 Hydrodynamics - Concepts and Experiments


[26] Gurel M., Ekdal A., Erturk A., Tanik A. Efforts towards setting eutrophication assess‐ ment criteria for coastal marine ecosystems. International Journal of Environment and Pollution 2005: 23(3), 325-335.

[14] Moncheva S., Gotsis-Skretas O., Pagou K., Krastev A. Phytoplankton blooms in Black Sea and Mediterranean coastal ecosystems subjected to anthropogenic eutrophica‐ tion: Similarities and differences. Estuarine, Coastal and Shelf Science 2001: 53(3),

[15] Herrera-Silveira J. A., Morales-Ojeda S. M. Evaluation of the health status of a coastal ecosystem in southeast Mexico: Assessment of water quality, phytoplankton and

[16] Taebi S., Etemad-Shahidi A., Faradi G. A. Examination of three eutrophication indi‐ ces to characterize water quality in the north east of Persian Gulf. J. Coastal Res. 2005:

[17] Zoriasatein N., Jalili S., Poor F. Evaluation of ecological quality status with the tro‐ phic index (TRIX) values in coastal area of Arvand, northeastern of Persian Gulf,

[18] Primpas I., Karydis M. Scaling the trophic index (TRIX) in oligotrophic marine envi‐ ronments. Environmental Monitoring and Assessment 2011: 178(1-4), 257-269.

[19] Toms J. D., Lesperance M. L. Piecewise regression: A tool for identifying ecological

[20] Nikolaidis G., Patoucheas D. P., Moschandreou K. Estimating breakpoints of Chl-a in relation with nutrients from Thermaikos Gulf (Greece) using piecewise linear regres‐

[21] Gaston A. J., Smith P. A., Provencher J. F. Discontinuous change in ice covers in Hud‐ son Bay in the 1990s and some consequences for marine birds and their prey. ICES

[22] Prasad A. K., Chai L., Singh R. P., Kafatos M. Crop yield estimation model for Iowa using remote sensing and surface parameters. International Journal of Applied Earth

[23] Cavanaugh K. C., Kellner J. R., Forde A. J., Gruner D. S., Parker J. D., Rodriguez W., Feller I. C. Poleward expansion of mangroves is a threshold response to decreased frequency of extreme cold events. Proceedings of the National Academy of Sciences

[24] Morley S. A., Lai C. H., Clarke A., Tan K. S., Thorne M. A., Peck L. S. Limpet feeding rate and the consistency of physiological response to temperature. Journal of Compa‐

[25] Sany S. B. T., Hashim R., Rezayi M., Salleh A., Safari O. A review of strategies to monitor water and sediment quality for a sustainability assessment of marine envi‐

ronment. Environmental Science and Pollution Research 2014: 21(2), 813-833.

sion. Fresenius Environmental Bulletin 2006: 15(9 B), 1189-1192.

Journal of Marine Science: Journal du Conseil 2012: fss 040.

Observation and Geoinformation 2006: 8(1), 26-33.

of the USA 2014: 111(2), 723-727.

rative Physiology B 2014: 1-8.

submerged aquatic vegetation. Marine Pollution Bulletin 2009: 59(1), 72-86.

281-295.

48 Hydrodynamics - Concepts and Experiments

42, 405-411.

Iran. World 2013: 5(3), 257-262.

thresholds. Ecology 2003: 84(8), 2034-2041.


[54] Yamamoto T., Oh S. J., Kataoka Y. Growth and uptake kinetics for nitrate, ammoni‐ um and phosphate by the toxic dinoflagellate *Gymnodinium catenatum* isolated from Hiroshima Bay, Japan. Fisheries Science 2004: 70(1), 108–115.

[40] Bissinger J. E., Montagnes D. J. S., Sharples J., Atkinson D. Predicting marine phyto‐ plankton maximum growth rates from temperature: Improving on the Eppley curve

[41] Steele J. H. Environmental control of photosynthesis in the sea. Limn. Oceanogr.

[42] Anonymous. Observations meteorologiques de Thessaloniki 2000-2001. Ann. de'l

[43] Putt M., Rivkin R., Prezelin B. B. Effects of altered photic regimes on diel patterns of species-specific photosynthesis - I. Comparison of polar and temperate phytoplank‐

[44] Andersen P., Nival P. Modelling of phytoplankton population dynamics in an en‐

[45] Landry M. R., Hassett R. P. Estimating the grazing impact of marine micro-zooplank‐

[46] Landry M. R., Peterson W. K., Lorenzent C. J. Zooplankton grazing, phytoplankton growth, and export flux: Inferences from chlorophyll tracer methods. ICES J. Mar.

[47] Peperzak L., Duin R. N. M., Colijn F., Gieskes W. W. C. Growth and mortality of flag‐ ellates and non-flagellate cells of *Phaeocystis globosa* (Prymnesiophyceae). Journal of

[48] Garcés E., Vila M., Masó M., Sampedro N., Giacobbe M. G., Penna A. Taxon-specific analysis of growth and mortality rates of harmful dinoflagellates during bloom con‐

[49] Huret M., Gohin F., Delmas D., Lunven M., Garçon V. G. Use of SeaWiFS data for light availability and parameter estimation of a phytoplankton production model of

[50] Eppley R. W., Rogers J. N., McCarthy J. J. Half saturation constants for uptake of ni‐ trate and ammonium by various phytoplankton. Limnol. Oceanogr. 1969: 14, 912-920.

[51] Pasciak W. J., Gavis J. Transport limitation of nutrient uptake in phytoplankton. Lim‐

[52] Taft J. L., Taylor W. R., McCarthy J. J. Uptake and release of phosphorus by phyto‐ plankton in the Chesapeake Bay estuary, USA. Marine Biology 1975: 33, 21-32. [53] Yamamoto T., Tarutani K. Growth and phosphate uptake kinetics of the toxic dino‐ flagellate *Alexandrium tamarense* from Hiroshima Bay in the Seto Inland Sea, Japan.

closed water column. J. Mar. Biol. Ass. UK 1989: 69, 625-646.

ditions. Marine Ecology Progress Series 2005: 301, 67-79.

the Bay of Biscay. Journal of Marine Systems 2007: 65, 509–531.

using quantile regression. Limnology and Oceanography 2008: 53(2), 487-493.

1962: 7, 137-150.

50 Hydrodynamics - Concepts and Experiments

Inst. Meteir. et Climatol. 2001: 69.

ton. Marine Biology 1988: 97, 435-443.

ton. Marine Biology 1982: 67, 283-288.

Plankton Research 2000: 22(1), 107-120.

nol. Oceanogr. 1974: 19, 881-885.

Phycological Research 1999: 47(1), 27-32.

Sci. 1995: 52, 337-345.


## **Alternative Uses of Cavitating Jets**

José Gilberto Dalfré Filho, Maiara Pereira Assis and Ana Inés Borri Genovez

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/59136

#### **1. Introduction**

[68] Rosenberg R. Benthic macrofaunal dynamics, production, and dispersion in an oxy‐ gen-deficient estuary of west Sweden. Journal of Experimental Marine Biology and

[69] Zarkanellas A. J. The effects of pollution-induced oxygen deficiency on the benthos in Elefsis Bay, Greece. Marine Environmental Research, 1979: 2(3), 191-207.

[70] Bogdanos C., Satsmadjis J. The macrozoobenthos of an Aegean embayment. Thalas‐

[71] Lie U., Kisker D. S. Species composition and structure of benthic infauna communi‐ ties off the coast of Washington. Journal of the Fisheries Board of Canada 1970:

[72] Diapoulis A., Bogdanos C. Preliminary study of soft substrate macrozoobenthos and marine flora in the Bay of Gera (Lesvos Island, Greece). Thalassographica, 1983: 6,

Ecology 1977: 26(2), 107-133.

sographica, 1983: 6, 77-105.

27(12), 2273-2285.

52 Hydrodynamics - Concepts and Experiments

127-139.

Cavitation is a phase transformation occurring in a fluid system under certain conditions. If the dynamic alteration of the absolute static pressure reaches or drops below the vapour pressure of the liquid, vapour bubbles are formed inside the liquid and can collapse as they move to a high-pressure region [1]. When the vapour bubbles collapse, shock waves are produced that propagate at the speed of sound through the liquid [2]. So, cavitation can produce undesirable effects such as noise, vibration, pressure fluctuation, erosion and efficiency loss in a hydraulic system.

The pressure difference due to the dynamic effect of the fluid motion is proportional to the square of the relative velocity. This can be written in the form of the usual pressure coefficient (Cp), presented in the Equation 1 [3]:

$$C\_p = \frac{(p - p\_o)}{\rho \, v\_o^2 / 2} \tag{1}$$

where ρ is the density of the liquid; vo and po are the velocity and pressure of an undisturbed liquid, respectively, and (p-po) is the pressure differential due to dynamic effects of fluid motion.

If pressure p reaches a minimum, p=pmin (see Figure 1), the pressure coefficient will be minimum, and a set of conditions can be created so that pmin drops to a value where cavitation can begin. This can be accomplished by raising the relative velocity vo for a fixed value of the pressure po or lowering po with vo remaining constant. If surface tension is ignored, the pressure pmin will be the pressure inside the cavity. This pressure will be the bubble pressure.

© 2015 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and eproduction in any medium, provided the original work is properly cited.

If we consider that cavitation will occur when the normal stresses at a point in the liquid are reduced to zero, then the bubble pressure will assume the value of the vapour pressure, pv. Then, the cavitation index (σ) is the theoretical value of the negative pressure coefficient, Equation 2:

$$\sigma = \frac{p\_o - p\_v}{\rho \, v\_o^2 / 2} = - (\mathbf{C}\_p)\_{\text{min}} \tag{2}$$

The cavitation index can designate the probability of a system to cavitate and establish different intensity levels of cavitation in this system.

**Figure 1.** Pressure distribution in a submerged body [4].

As the pressure is increased, the bubble diameter decreases from the original size reaching a minimum size. The process takes place until the bubble diameter becomes microscopic. If the bubble collapses, then shock waves form with celerity equal to the speed of sound in water. If a boundary is close to where the bubble imploded, it will deform into a microjet. The velocity of this microjet is high and the shockwave produces high pressure responsible for the cavita‐ tion damage of a surface [5]. The collapse of the cavities formed by the cavitating jet generates high-pressure waves, estimated to be approximately 69.00 MPa [6] and high-speed microjets (above 100 m/s), all of which have a significant amount of destructive power [7].

Flow conditions leading to the onset of cavitation are generally conservative in predicting damage. The severity of the damage may be related to both intensity of cavitation and exposure time [8].

In hydraulic structures, as the high velocity flow passes over the many irregularities that exist in the concrete surface, cavitation can commence and consequently, damages may occur. Therefore, the material properties of the surface have to be improved to provide adequate resistance either during the construction phase or when the substitution of the eroded concrete is required.

Particularly of interest, in developing countries, hydropower will continue to play a significant role in supplying energy. Therefore, as energy demands are continuously increasing, new hydropower plants (dams and appurtenances) are being constructed. Nevertheless, the operation and maintenance of these structures, which include the spillways and tunnels under high velocity flows, are a great issue. Therefore, the material properties have to be improved to provide adequate resistance, either during the construction phase or substitution phase. For example, the substitution of a great area of eroded concrete was required in the case of Porto Colombia hydropower spillway and dissipation basin in the Grande River, Brazil [7].

An additional complexity is observed in hydraulic structures due to the simultaneous effect of cavitation and high impact of the flow [9,10,4]. Evaluation of the concrete resistance erosion in hydraulic structures is essential to guarantee adequate operation. It was suggested by [11] defining a methodology to appropriately test the materials, when submitted to cavitation, to be used in hydraulic structures.

In this chapter the authors discuss the use of high velocity cavitating jets to determine the erosion in high performance concretes for hydraulic structures. Moreover, two alternatives uses for cavitating jets are presented: 1) the inactivation of *Escherichia coli,* and, 2) the decom‐ position of persistent compounds in water.

#### **2. Use of the cavitating jet to test erosion in surfaces**

#### **2.1. Introduction**

If we consider that cavitation will occur when the normal stresses at a point in the liquid are reduced to zero, then the bubble pressure will assume the value of the vapour pressure, pv. Then, the cavitation index (σ) is the theoretical value of the negative pressure coefficient,

<sup>2</sup> min ( ) <sup>2</sup>

The cavitation index can designate the probability of a system to cavitate and establish different

As the pressure is increased, the bubble diameter decreases from the original size reaching a minimum size. The process takes place until the bubble diameter becomes microscopic. If the bubble collapses, then shock waves form with celerity equal to the speed of sound in water. If a boundary is close to where the bubble imploded, it will deform into a microjet. The velocity of this microjet is high and the shockwave produces high pressure responsible for the cavita‐ tion damage of a surface [5]. The collapse of the cavities formed by the cavitating jet generates high-pressure waves, estimated to be approximately 69.00 MPa [6] and high-speed microjets

Flow conditions leading to the onset of cavitation are generally conservative in predicting damage. The severity of the damage may be related to both intensity of cavitation and exposure

In hydraulic structures, as the high velocity flow passes over the many irregularities that exist in the concrete surface, cavitation can commence and consequently, damages may occur. Therefore, the material properties of the surface have to be improved to provide adequate resistance either during the construction phase or when the substitution of the eroded concrete

(above 100 m/s), all of which have a significant amount of destructive power [7].

*p*


*o v*

s

intensity levels of cavitation in this system.

**Figure 1.** Pressure distribution in a submerged body [4].

*o p p <sup>C</sup> <sup>v</sup>*

r

Equation 2:

54 Hydrodynamics - Concepts and Experiments

time [8].

is required.

Evaluation of the concrete erosion resistance in hydraulic structures is essential to guarantee adequate operation, and tests should be performed under the same operation conditions. Various techniques can be used to induce cavitation, such as, ultrasonic methods, hydrody‐ namic methods, and high-speed/high-pressure homogenization. In the hydrodynamic cavitation, pressure variations are produced using the geometry of the system, while in the acoustic cavitation pressure variations are effected using sound waves.

Concrete resistance to cavitation damage was investigated by [12-14] using a Venturi device. A chamber was used by [15] to test concrete samples subjected to short-duration cavitation. An alternative test for erosion evaluation uses water cavitating jet technology [16,17] to achieve short test time. A flow cavitation chamber was used by [18-20] to test concrete and different rock samples subjected to short-duration cavitation. Water cavitating jet technology was firstly used to clean surfaces. For example, [21] compared the efficiency of conventional noncavitating jets to cavitating jets.

Regardless several researchers have worked with the cavitation jet apparatus, to test concretes submitted to the erosive effect of cavitation, the appropriate nozzle specifications and sample dimensions must be experimentally determined.

#### **2.2. The cavitating jet apparatus**

The cavitating jet apparatus uses a nozzle specially designed to produce cavitation, combining high-velocity flows and cavitation (with an appropriate cavitation index). Because the collapse of bubbles is concentrated over a microscopic area, localized stresses are produced, providing the cavitating jet with a great advantage over steady non-cavitating jet operating at the same pressures and flow rates [22]. Figure 2 shows a schematic of the formation of a cavitating jet. As the high velocity flow leaves the orifice nozzle, eddies are formed between the high velocity layer of the jet and the surrounding liquid. If a nucleus in the water is captured in one of these eddies and the pressure inside the eddy drops to vapour pressure, then the nucleus will begin to grow. However, if the pressure remains near the vapour pressure long enough for the nuclei to reach the critical diameter, then it begins to grow almost by vaporization. As long as the size of the vapour cavity increases, the strength of the eddy is destroyed and the rotational velocity decreases. Then, the surrounding bubble pressure is no longer the vapour pressure. Inside the cavity the pressure remains at vapour pressure and cavity surroundings at hydro‐ static pressure. Consequently, the cavity becomes unstable and it collapses inward.

**Figure 2.** Schematic flow of a submerged jet [6].

Now, consider in Figure 3 a spherical vapour bubble of initial radius b and internal pressure Po. Later on, the bubble acquires the radius re at pressure Pe. If a small amount of air is in the liquid, such that the cavities are filled exclusively by vapour, the bubble growth and collapse is intense and causes severe damages to the vicinity [6].

The pressure magnitude generated at the solid boundary is expected to be a main factor in estimating the erosion efficiency of a cavitating jet impinging on a solid boundary as shown in Figure 3.

For a steady water jet (vj ),the generated pressure magnitude (Ps) can be estimated from the stagnation pressure, as in Equation 3:

$$P\_s = \bigvee\_2 \rho \mathbf{v}\_j^2 \tag{3}$$

A moderate level pressure would be inappropriate for testing erosion rates. Depending on the standoff distance, water droplets formation can occur and the impact pressure can be larger

**Figure 3.** Cavitating jet formation and bubble collapse (Adapted from [23])

than Ps. It was indicated that the transient pressure (Pa) between the cavity implosion and the solid boundary could be approximated by the "water hammer" equation, as in Equation 4:

$$P\_{\mathbf{a}} = \rho \mathbf{c} \mathbf{v}\_{j} \tag{4}$$

with c as the sound velocity in water.

**2.2. The cavitating jet apparatus**

56 Hydrodynamics - Concepts and Experiments

**Figure 2.** Schematic flow of a submerged jet [6].

stagnation pressure, as in Equation 3:

in Figure 3.

For a steady water jet (vj

is intense and causes severe damages to the vicinity [6].

The cavitating jet apparatus uses a nozzle specially designed to produce cavitation, combining high-velocity flows and cavitation (with an appropriate cavitation index). Because the collapse of bubbles is concentrated over a microscopic area, localized stresses are produced, providing the cavitating jet with a great advantage over steady non-cavitating jet operating at the same pressures and flow rates [22]. Figure 2 shows a schematic of the formation of a cavitating jet. As the high velocity flow leaves the orifice nozzle, eddies are formed between the high velocity layer of the jet and the surrounding liquid. If a nucleus in the water is captured in one of these eddies and the pressure inside the eddy drops to vapour pressure, then the nucleus will begin to grow. However, if the pressure remains near the vapour pressure long enough for the nuclei to reach the critical diameter, then it begins to grow almost by vaporization. As long as the size of the vapour cavity increases, the strength of the eddy is destroyed and the rotational velocity decreases. Then, the surrounding bubble pressure is no longer the vapour pressure. Inside the cavity the pressure remains at vapour pressure and cavity surroundings at hydro‐

static pressure. Consequently, the cavity becomes unstable and it collapses inward.

Now, consider in Figure 3 a spherical vapour bubble of initial radius b and internal pressure Po. Later on, the bubble acquires the radius re at pressure Pe. If a small amount of air is in the liquid, such that the cavities are filled exclusively by vapour, the bubble growth and collapse

The pressure magnitude generated at the solid boundary is expected to be a main factor in estimating the erosion efficiency of a cavitating jet impinging on a solid boundary as shown

> 1 <sup>2</sup> <sup>2</sup> *s j <sup>P</sup>* <sup>=</sup> rn

A moderate level pressure would be inappropriate for testing erosion rates. Depending on the standoff distance, water droplets formation can occur and the impact pressure can be larger

),the generated pressure magnitude (Ps) can be estimated from the

(3)

However, the pressure magnitude generated with the impact of a cavitating jet (Pb) is different from that generated by the impact of distinct droplets. According to [21] the maximum pressure Pb develops if a bubble collapses in an incompressible flow, under isothermal conditions, resulting in equation (5):

$$P\_b = \frac{P\_S}{6\epsilon\,35}e^{\left(\frac{\lambda}{\mu}\right)}\tag{5}$$

where α relates the gas pressure inside the bubble in the beginning of collapse under the pressure Ps. Therefore, the cavitating jet erosion efficiency R\* is related to Pb and Pa in Equation (6):

$$R^\star = \frac{P\_b}{P\_a} \tag{6}$$

Equation (5), however, does not take into account the geometry and finishing of the nozzle. It was verified by [10] that the finish and the geometry of the nozzle influence the erosion performance of the equipment. The author accomplished a series of tests with a cavitating jet apparatus, verifying the influence of these two variables in aluminium and concrete samples. Therefore, in Equation (7) a geometric efficiency factor η was introduced in the cavitation erosion efficiency [7] to calculate the total efficiency η :

$$
\varepsilon = \eta R \,\,\mathrm{\*}\tag{7}
$$

#### **2.3. Experimental setup**

To simulate the cavitation phenomenon and to evaluate the erosion in samples, a cavitating jet apparatus was constructed [10] in the Laboratory of Hydraulics and Fluid Mechanics at the State University Campinas, Brazil. Figure 5 and 6 shows a schematic representation and a photo of the test facility. A high-pressure displacement pump (kept at pressure 12.00 MPa) conducts water to the facility and to the metallic chamber provided with two windows for visualization of the tests. A pressure-regulating valve guaranteed safe operation. At the end of the high-pressure pipe and inside of the inactivation chamber, the nozzle with a hole (diameter of 2.00 x 10-3 m) was positioned at the end of the high-pressure pipe.

The nozzle specifications must be established according to its use and they are determined experimentally. The researcher [10] used three nozzle geometries (20° conical, 132° conical and circular) with both shaped and rounded edges (Figure 4).

**Figure 4.** Nozzle geometries (a) conical 20°, (b) conical 132°, (c) circular [7].

The facility is supplied with no re-circulated water. The samples are placed in the chamber subjected to the cavitating jet, that is, the jet is impinging onto the sample surface. The chamber is filled with water in order to allow the occurrence of the bubbles implosion. The damaged area was measured, the pits number were counted and the samples were photographed every 60 seconds until 300 seconds of time test. Then, the test continued more 900 uninterrupted seconds, achieving a total test time of 1200 seconds.

Apparatus test parameters can be adjusted to obtain optimal efficiency. In this experiment, pressure and velocity of the system were adjusted to reach a cavitation index of 0.14, which is considered to cause damage in hydraulic structures [11]. The three different geometry nozzles were used to evaluate the efficiency of cavitating jets.

**Figure 5.** Schematic of experimental setup for erosion tests.

a \* *<sup>b</sup> <sup>P</sup> <sup>R</sup>*

Equation (5), however, does not take into account the geometry and finishing of the nozzle. It was verified by [10] that the finish and the geometry of the nozzle influence the erosion performance of the equipment. The author accomplished a series of tests with a cavitating jet apparatus, verifying the influence of these two variables in aluminium and concrete samples. Therefore, in Equation (7) a geometric efficiency factor η was introduced in the cavitation

To simulate the cavitation phenomenon and to evaluate the erosion in samples, a cavitating jet apparatus was constructed [10] in the Laboratory of Hydraulics and Fluid Mechanics at the State University Campinas, Brazil. Figure 5 and 6 shows a schematic representation and a photo of the test facility. A high-pressure displacement pump (kept at pressure 12.00 MPa) conducts water to the facility and to the metallic chamber provided with two windows for visualization of the tests. A pressure-regulating valve guaranteed safe operation. At the end of the high-pressure pipe and inside of the inactivation chamber, the nozzle with a hole

The nozzle specifications must be established according to its use and they are determined experimentally. The researcher [10] used three nozzle geometries (20° conical, 132° conical and

e h

(diameter of 2.00 x 10-3 m) was positioned at the end of the high-pressure pipe.

erosion efficiency [7] to calculate the total efficiency η :

circular) with both shaped and rounded edges (Figure 4).

**Figure 4.** Nozzle geometries (a) conical 20°, (b) conical 132°, (c) circular [7].

**2.3. Experimental setup**

58 Hydrodynamics - Concepts and Experiments

*<sup>P</sup>* <sup>=</sup> (6)

= *R* \* (7)

Concrete is a heterogeneous material, which poses difficulties to test and evaluate optimal conditions based on comparison test. Thus, aluminum samples were used to establish the nozzle efficiency because the metal is homogenous and allows a comparison among tests in order to evaluate the differing optimal erosion conditions, and damages. Cavitation erosion rates can be analyzed in terms of number of pits by time in order to obtain quantitative and qualitative information on the erosion intensity variation. A binocular magnifier 160x was used to count pits in the surface of the aluminum samples through the tests

First, the apparatus tests were run with the different nozzles presented in Figure 4. The typical value of α=1/4 used by [21] was adopted and maintained for comparison purposes. The sound velocity c in water at 20°C was assumed to be 1482 m/s and water density as 998.2 kg/m3 . Given the value of the jet velocity vj, calculations were done for all the data and the efficiencies η were calculated. It was presented by [7] the efficiency results for the nozzles (Table1).

**Figure 6.** Photo of experimental setup in the Laboratory.


**Table 1.** Total erosion efficiency (ε) [7]

According to Table 1, better results were observed for sharp-edged nozzles than for chamferededge ones. The best total efficiency was obtained for the 132° conical nozzle (ε=0.43). Also, the time intervals required to accomplish the same erosion rate were observed to be 300 seconds for 132° conical sharp-edged nozzle, and 1200 seconds for the chamfered-edge nozzle. Therefore, not only the entrance pressure, but also nozzle geometry and finish are important for optimizing the apparatus efficiency. The authors pointed out that obtuse-angle nozzles with sharp-edge are preferable to perform cavitation erosion tests.

Using cavitating jet apparatus (σ=0.14) and the 132° conical sharp-edged nozzle, three special concrete samples were tested to obtain the cavitation erosion rate. The procedure of making samples begins by separating, cleaning and stocking the aggregates in the laboratory the day before. The next day, samples were molded in agreement with the Brazilian Standard NBR 5738. Table 2 presents samples compositions and compressive resistance obtained.


**Table 2.** Characteristics of the concrete samples.

Figure 7 shows the erosion comparison between concrete samples in terms of volume erosion in time.

**Figure 7.** Concrete erosion over time.

**Nozzle Sharp-edge (ε) Chamfered- edge (ε)**

According to Table 1, better results were observed for sharp-edged nozzles than for chamferededge ones. The best total efficiency was obtained for the 132° conical nozzle (ε=0.43). Also, the time intervals required to accomplish the same erosion rate were observed to be 300 seconds for 132° conical sharp-edged nozzle, and 1200 seconds for the chamfered-edge nozzle. Therefore, not only the entrance pressure, but also nozzle geometry and finish are important for optimizing the apparatus efficiency. The authors pointed out that obtuse-angle nozzles

Using cavitating jet apparatus (σ=0.14) and the 132° conical sharp-edged nozzle, three special concrete samples were tested to obtain the cavitation erosion rate. The procedure of making samples begins by separating, cleaning and stocking the aggregates in the laboratory the day before. The next day, samples were molded in agreement with the Brazilian Standard NBR

5738. Table 2 presents samples compositions and compressive resistance obtained.

132°Conical 0.43 0.17 Circular 0.33 0.12 20°Conical 0.16 0.12

with sharp-edge are preferable to perform cavitation erosion tests.

**Table 1.** Total erosion efficiency (ε) [7]

60 Hydrodynamics - Concepts and Experiments

**Figure 6.** Photo of experimental setup in the Laboratory.

The best results were obtained for sample C which contains hard aggregate, superior axial compressive resistance and the addition of silica. However sample A (fck 55.00 MPa) presented better results than sample B (fck 63.00 MPa), despite of the greater compressive resistance, showing that an adequate concrete (resistant to cavitation erosion) is a combination of several factors such as aggregate type, size and shape, additions to cement and water cement relation. Figure 8 shows the evolution of erosion in the B concrete sample. Additional tests were performed with unsubmerged and submerged samples, showing that the formation of a highvelocity cavitating jet provoked high erosion rates in the sample as expected, according [7]. In [24] the authors compared the erosion generated by the cavitating jet imping directly in to and parallel to a concrete sample. In this last case, the sample is not subjected to the effects of the impact force of the jet. Comparing the results, a higher erosion rate was observed when the samples are positioned directly over to the cavitating jet (Figure 9).

**Figure 8.** Erosion evolution of concrete sample B over time.

**Figure 9.** Concrete sample placed parallel to the cavitating jet [24].

Fast, efficient and economic testing is needed. The use of the cavitating jet allowed significant reduction in testing times, [10] specially when compared to the Venturi device. In [12] the author took 30 hours to compare the cavitating wear among concrete samples. As concrete is a heterogeneous material, samples are necessary to be larger than a few millimeters, as used by [25] to test metal samples. [10] used samples 20 cm in diameter and 5 cm in height to test concrete samples. The cavitating jet apparatus used in this experiment is compact, low cost and has short test times. It can also be used to measure erosion for testing concrete use in hydraulic structures. The cavitating jet apparatus creates a force larger than the one generated by a simple jet of high pressure, and thus can simulate the combined effect of high-speed flows and cavitation normally experienced in hydraulic structures

## **3. Use of the cavitating jet to inactivate bacteria and to decompose persistent compounds in water**

#### **3.1. Introduction**

**Figure 8.** Erosion evolution of concrete sample B over time.

62 Hydrodynamics - Concepts and Experiments

Water quality has deteriorated over the years due to industrial and agricultural pollution as well as an increase in the domestic sewage generated by a rapidly growing population. Pathogenic microorganisms in water pose a serious health and security threat to drinking water supply systems. The quality of water for human consumption can be improved by controlling pollution and by increasing the efficiency of inactivation techniques, which involves the destruction of pathogens present in water at a reasonable cost.

The use of chemicals, such as chlorine gas, sodium hypochlorite, calcium hypochlorite, hydrated ammonia, ammonium hydroxide, ammonium sulphate, and ozone, is common during the water disinfection process in conventional water treatment systems [26]. Many problems arise when using chemical methods for disinfection, including the high maintenance demands of the associated facilities (corrosion, incrustation), the formation of toxic byproducts (chlorine addition may generate byproducts such as trihalomethanes) and environmental concerns (chemical effluents released into rivers compromise aquatic life) [27]. Therefore, the development of alternative techniques for inactivating pathogens in water is desirable. Inactivation based on the cavitation phenomenon appears to be a promising alternative or supplement to existing techniques [28].

Research using the phenomenon of cavitation to inactivate microorganisms has been per‐ formed in the United States, Russia, India, China, Japan, UK, South Africa, France, Mexico, among others, using different methods. The hydrodynamic cavitation method is the most studied and disseminated. This is because the best results are achieved. However, there is currently no standardized method to carry out inactivation with this technique with respect to the pressures and speeds testing.

The performance shown by the studies on the wear of concrete indicated that the cavitating jet apparatus could be adapted to test the inactivation of microorganisms and the decompo‐ sition of persistent compounds in water. However, tests are necessary to establish appropriate pressures and velocities for each study. The inactivation of microorganisms was performed with *Escherichia coli* that are microorganisms commonly found in the intestines of humans and other warm-blooded animals and indicate the presence of water contamination. The decom‐ position of persistent compounds was performed with methylene blue.

#### **3.2. Experimental tests for inactivation and chemical decomposition**

As illustrated in Figure 2 the cavitating jet apparatus with 132° conical nozzle was adapted to conduct the experiments to inactivate *Escherichia coli* and to decompound the persistent methylene blue compound. The chamber used was 700.00 x 10-3 m high with a diameter of 300.00 x 10-3 m. Part of the container was filled with 40.00 x 10-3 m3 of contaminated water to be studied in a closed circuit. A low-cost refrigeration system was adapted to the apparatus. A 12.70 x 10-3 m diameter and 19,500.00 x 10-3 m long copper coil was placed inside of the inactivation chamber. This coil was connected to a reservoir with a capacity of 750.00 x 10-3 m3 that was filled with clean water. The temperature of the tests was held constant at 33 °C, thus assuring that the inactivation was exclusively due to the cavitation phenomenon instead of an increase in the water temperature.

All of the microbiological, physical, and chemical procedures of this study follow the proce‐ dures described in the Standard Methods for Examination of Water and Wastewater Manual [19,20]. Non-pathogenic *Escherichia coli* (ATCC 25922) were used for the microbiological tests. Viable cell counting was performed using the Colilert Method®. The presented results of the microbiological tests are the average of the results obtained in the dilutions.

The high-pressure displacement pump recirculated the water to be treated and was kept at a pressure of 4.00, 8.00, 10.00 and 12.00 MPa. First, the inactivation device was turned on without pressure to circulate the residue being treated. Samples were then collected at the initial time (T0) to obtain a control sample. Next, a vial of the sample was removed to measure the number of viable *Escherichia coli* cells*.* In addition to the control time point (T0), samples were collected every 900 seconds. All of the tests were repeated in triplicate, however, the results represent an average of the individual results.

Figure 10 shows that the higher the pressure, the higher the inactivation. At time test 1800 seconds, the following inactivation rates are observed: less than 37.50% for pressure test 4.00 MPa, 98.30% for test pressure 8.00 MPa, 99.96% for pressure test 10.00 MPa, and 100% for pressure 12.00 MPa,. Although the inactivation percentages are close to pressures 8.00, 10.00 and 12.00 MPa, the highest inactivation rate was achieved for test pressure 12.00 MPa.

**Figure 10.** Inactivation rate and pressure.

controlling pollution and by increasing the efficiency of inactivation techniques, which

The use of chemicals, such as chlorine gas, sodium hypochlorite, calcium hypochlorite, hydrated ammonia, ammonium hydroxide, ammonium sulphate, and ozone, is common during the water disinfection process in conventional water treatment systems [26]. Many problems arise when using chemical methods for disinfection, including the high maintenance demands of the associated facilities (corrosion, incrustation), the formation of toxic byproducts (chlorine addition may generate byproducts such as trihalomethanes) and environmental concerns (chemical effluents released into rivers compromise aquatic life) [27]. Therefore, the development of alternative techniques for inactivating pathogens in water is desirable. Inactivation based on the cavitation phenomenon appears to be a promising alternative or

Research using the phenomenon of cavitation to inactivate microorganisms has been per‐ formed in the United States, Russia, India, China, Japan, UK, South Africa, France, Mexico, among others, using different methods. The hydrodynamic cavitation method is the most studied and disseminated. This is because the best results are achieved. However, there is currently no standardized method to carry out inactivation with this technique with respect

The performance shown by the studies on the wear of concrete indicated that the cavitating jet apparatus could be adapted to test the inactivation of microorganisms and the decompo‐ sition of persistent compounds in water. However, tests are necessary to establish appropriate pressures and velocities for each study. The inactivation of microorganisms was performed with *Escherichia coli* that are microorganisms commonly found in the intestines of humans and other warm-blooded animals and indicate the presence of water contamination. The decom‐

As illustrated in Figure 2 the cavitating jet apparatus with 132° conical nozzle was adapted to conduct the experiments to inactivate *Escherichia coli* and to decompound the persistent methylene blue compound. The chamber used was 700.00 x 10-3 m high with a diameter of 300.00 x 10-3 m. Part of the container was filled with 40.00 x 10-3 m3 of contaminated water to be studied in a closed circuit. A low-cost refrigeration system was adapted to the apparatus. A 12.70 x 10-3 m diameter and 19,500.00 x 10-3 m long copper coil was placed inside of the inactivation chamber. This coil was connected to a reservoir with a capacity of 750.00 x 10-3 m3 that was filled with clean water. The temperature of the tests was held constant at 33 °C, thus assuring that the inactivation was exclusively due to the cavitation phenomenon instead of an

All of the microbiological, physical, and chemical procedures of this study follow the proce‐ dures described in the Standard Methods for Examination of Water and Wastewater Manual [19,20]. Non-pathogenic *Escherichia coli* (ATCC 25922) were used for the microbiological tests.

position of persistent compounds was performed with methylene blue.

**3.2. Experimental tests for inactivation and chemical decomposition**

involves the destruction of pathogens present in water at a reasonable cost.

supplement to existing techniques [28].

64 Hydrodynamics - Concepts and Experiments

to the pressures and speeds testing.

increase in the water temperature.

However, the optimum inactivation rate will be achieved by relating the inactivation rate to the energy consumption. The Energy Efficiency (*EE*) of the apparatus is calculated with Equation 8 and the experiment results are shown in Table 3.

$$EE = \frac{Ci - Cf \cdot \forall}{P \cdot T} \tag{8}$$

being:

Cf: Final Concentration of *Escherichia coli* (CFU/mL);

Ci: Initial Concentration of *Escherichia coli* (CFU/mL);

∀: Volume (mL);

P: Power (W);

T: Time (seconds).


**Table 3.** Energy Efficiency to inactivate *Escherichia Coli* at different pressure.

In order to compare the tests, a time test of 1800 seconds was standardized. Table 3 indicates that for the conical 132° nozzle, the most efficient operation pressure is 10.00 MPa, followed by pressure 4.00, 12.00 and 8.00 MPa.

A complementary study was conducted with natural waters retrieved from a lake at the University's campus [28]. Given past studies by [29,30], high efficiency for the first 900 seconds of the test with an inactivation rate of 90% was obtained.

The substance methylene blue was used for testing degradation of persistent compounds. The results for each time and pressure test were analyzed using the spectrophotometer (Hach DR / 4000U Spectrophotometer) scan from 800 nm to 200 nm. To analyze the degradation of persistent compounds it is necessary to observe the decrease in absorbance at a given wave‐ length band over time. The absorbance is proportional to the concentration of the substance based on the law of Lambert-Bee. Changes in the molecule of methylene blue are then analyzed for two wavelength bands: the ranges from 200nm to 400nm, ultraviolet (UV), and the range from 400nm to 800nm, visible to the human eye.

The result of the degradation is shown in Figures 11 and 12 at a pressure of 12.00 MPa and time 1800 seconds. Figure 11 shows no significant reduction in absorbance. Figure 12 indicates a reduction of 10.5% in the final absorbance at the peak wavelength range of 664nm.

**Figure 11.** Decomposition of methylene blue, pressure 12 MPa, UV wavelength.

being:

∀: Volume (mL);

66 Hydrodynamics - Concepts and Experiments

T: Time (seconds).

P: Power (W);

Cf: Final Concentration of *Escherichia coli* (CFU/mL); Ci: Initial Concentration of *Escherichia coli* (CFU/mL);

**Table 3.** Energy Efficiency to inactivate *Escherichia Coli* at different pressure.

of the test with an inactivation rate of 90% was obtained.

from 400nm to 800nm, visible to the human eye.

by pressure 4.00, 12.00 and 8.00 MPa.

**Parameters Test pressures**

**4.00 MPa 8.00 MPa 10.00 MPa 12.00 MPa** t (s) 1800.00 1800.00 1800.00 1800.00 V (mL) 40000.00 40000.00 40000.00 40000.00 i (A) 25.00 29.00 31.00 32.80 U (V) 220.00 220.00 220.00 220.00 P (W) 5500.00 6380.00 6820.00 7216.00 Ci (CFU/mL) 918425.00 141360.00 800013.00 256547.00 Cf (CFU /mL) 574000.00 2410.00 337.00 9.00 Ci-Cf (CFU /mL) 344425.00 138950.00 799676.00 256538.00 (Ci-Cf)/P(CFU /mL)/(W) 62.62 21.78 117.25 35.55 ЄЄ (CFU /J) 1391.62 483.98 2605.66 790.03

In order to compare the tests, a time test of 1800 seconds was standardized. Table 3 indicates that for the conical 132° nozzle, the most efficient operation pressure is 10.00 MPa, followed

A complementary study was conducted with natural waters retrieved from a lake at the University's campus [28]. Given past studies by [29,30], high efficiency for the first 900 seconds

The substance methylene blue was used for testing degradation of persistent compounds. The results for each time and pressure test were analyzed using the spectrophotometer (Hach DR / 4000U Spectrophotometer) scan from 800 nm to 200 nm. To analyze the degradation of persistent compounds it is necessary to observe the decrease in absorbance at a given wave‐ length band over time. The absorbance is proportional to the concentration of the substance based on the law of Lambert-Bee. Changes in the molecule of methylene blue are then analyzed for two wavelength bands: the ranges from 200nm to 400nm, ultraviolet (UV), and the range

**Figure 12.** Decomposition of methylene blue, pressure 12 MPa, visible wavelength.

Comparing the results of Figures 11 and 12, the greatest difference in absorbance occurs in the visible band. Pressures lower than 12.00 MPa did not show any significant decrease in the absorbance. At the end of the tests, the residue was placed in translucent plastic bottles and exposed to sunlight. It was observed that the methylene blue dye was then degraded only by the peroxidation generated by the cavitation treatment.

The synergy between the cavitation process and the oxidative processes is a current and relevant research [31-33]. Recent work by [34-36] has been focusing on the degradation of persistent pollutants, such as pharmaceuticals, dyes, industrial effluents, pesticides, in addition to research with microorganisms.

## **4. Final remarks**

Cavitation is known for its undesired effects produced in hydraulic systems, namely, the noise, vibrations, pressure fluctuations, erosion and efficiency loss. In hydraulic structures, cavitation can be destructive. Many times, the materials to do repairs in the structure are used without the necessary laboratory tests. These tests can be onerous, showing the importance of previous studies to check the applicability, either for the project phase or for repairs.

Despite the undesirable effects of cavitation, this phenomenon showed to be an interesting alternative to inactivate bacteria in water. The use of chemicals, such as chlorine, is commonly used for water disinfection in conventional water treatment systems. Cavitation can be used as a "clean" treatment, as it can reduce the quantity of chemicals added to water in the conventional water treatments.

The cavitating jet equipment simulates cavitation and its efficiency can be improved by using different nozzle configuration that forms the cavitating jet. This work showed that the most efficient nozzle geometry was obtained using a sharp-edge 132° conical nozzle for different applications, erosion tests, bacteria inactivation and decomposition of compounds.

The cavitating jet apparatus proved to be efficient, low cost, low power consumption when used to test material erosion and to inactivate *Escherichia coli*. Furthermore, it can be adapted to operate at larger scale treatments.

#### **Acknowledgements**

The authors would like to acknowledge FAPESP (São Paulo Research Foundation) for the Scientific Research Scholarship (Process: 2000/03611-0, 2000/03732-2, 2002/10348-0, 2009/53553-1 and 2011/16347-4) and Research Award Scheme (Process: 2009/54278-4).

#### **Author details**

Comparing the results of Figures 11 and 12, the greatest difference in absorbance occurs in the visible band. Pressures lower than 12.00 MPa did not show any significant decrease in the absorbance. At the end of the tests, the residue was placed in translucent plastic bottles and exposed to sunlight. It was observed that the methylene blue dye was then degraded only by

The synergy between the cavitation process and the oxidative processes is a current and relevant research [31-33]. Recent work by [34-36] has been focusing on the degradation of persistent pollutants, such as pharmaceuticals, dyes, industrial effluents, pesticides, in

Cavitation is known for its undesired effects produced in hydraulic systems, namely, the noise, vibrations, pressure fluctuations, erosion and efficiency loss. In hydraulic structures, cavitation can be destructive. Many times, the materials to do repairs in the structure are used without the necessary laboratory tests. These tests can be onerous, showing the importance of previous

Despite the undesirable effects of cavitation, this phenomenon showed to be an interesting alternative to inactivate bacteria in water. The use of chemicals, such as chlorine, is commonly used for water disinfection in conventional water treatment systems. Cavitation can be used as a "clean" treatment, as it can reduce the quantity of chemicals added to water in the

The cavitating jet equipment simulates cavitation and its efficiency can be improved by using different nozzle configuration that forms the cavitating jet. This work showed that the most efficient nozzle geometry was obtained using a sharp-edge 132° conical nozzle for different

The cavitating jet apparatus proved to be efficient, low cost, low power consumption when used to test material erosion and to inactivate *Escherichia coli*. Furthermore, it can be adapted

The authors would like to acknowledge FAPESP (São Paulo Research Foundation) for the Scientific Research Scholarship (Process: 2000/03611-0, 2000/03732-2, 2002/10348-0,

2009/53553-1 and 2011/16347-4) and Research Award Scheme (Process: 2009/54278-4).

applications, erosion tests, bacteria inactivation and decomposition of compounds.

studies to check the applicability, either for the project phase or for repairs.

the peroxidation generated by the cavitation treatment.

addition to research with microorganisms.

68 Hydrodynamics - Concepts and Experiments

**4. Final remarks**

conventional water treatments.

to operate at larger scale treatments.

**Acknowledgements**

José Gilberto Dalfré Filho\* , Maiara Pereira Assis and Ana Inés Borri Genovez

\*Address all correspondence to: dalfre@fec.unicamp.br

College of Civil Engineering, Architecture, and Urbanism, State University of Campinas, UNICAMP, Campinas, Brazil

#### **References**


[25] Lichtarowicz, A. Cavitating jet apparatus for cavitating erosion testing. Erosion: pre‐ vention and useful applications. American Society for Testing and Materials, Phila‐ delphia, STP 664, p. 530-549,1979.

[11] Semenkov, V.M. 2007. Discussion on the work "Evaluating damage of concrete hy‐

[12] Gikas, I. Cavitação-Efeitos sobre superfícies de resina de epóxi e concretos comuns e

[13] Billard, J.Y. and Fruman, D.H. (1992) "Etude expérimentale de i'influence d'un champ de pression fluctuant sur l'apparition de la cavitation dans un Venturi". La

[14] Pham, T. M.; Michel, J.M.; Lecoffre,Y. (1997) "Qualification et optimization d'un compteur dynamique à ogive centrale". La Houille Blanche, Paris, année 52, n.342, p.

[15] Momber, A.W. (2000) Short time cavitation erosion of concrete. WEAR, n.241, p.

[16] Conn et al. Some unusual applications for cavitating water jets. In: International Symposium on Jet Cutting Technology, 7, 1984, Ottawa. Jet Cutting Technology.

[17] Cheng, C.L.; Webster, C.T. e Wong, J.Y. Cavitation resistant coatings for concrete hy‐ draulic structures. ACI Materials Journal, Detroit, v.87, n.6, p.594-601, nov-dec1990.

[18] Momber, A.W. Deformation and fracture of rocks due to high-speed liquid impinge‐

[19] Momber, A.W. Aggregate liberation from concrete by flow cavitation. Int. J. Miner.

[20] Hu, X.G.; Momber A.W.; Yin, Y. Erosive wear of hydraulic concrete with low steel fiber content. Journal of Hydraulic Enmgineering. ASCE.. p 1331-1340. Dec. 2006. [21] Houlston, R., Vickers, G.W. (1978). Surface cleaning using water jet cavitation and droplet erosion. In: 4th Int. Sym. on Jet Cutting Technology. Cranfield: BHRA, The

[22] Kalumuck, K.M., Chahine, G.L., Hsiao, C.-T., Choi, J.-K. (2003). Remediation and dis‐ infection of water using jet generated cavitation. Proc. 5th Int. Sym. Cavitation, Osa‐

[23] Vickers, G.W., Houlston, R. Modelling the erosion efficiency of cleaning jets. Applied

[24] Dalfré Filho, J.G.; Genovez, A.I.B.; Luvizotto Junior (2006). Equipamento para Simu‐ lação da Erosão por Cavitação. Revista Brasileira de Recursos Hídricos, v.11, n.1, p.

Cranfield: BHRA The Fluid Engineering Centre, june,1984. p.1-12.

ment. International Journal of Fracture 130: 683–704, 2004.

especiais. Boletim técnico DAEE, São Paulo,4(1), p.89-121, jan.-abril 1981.

draulic structures caused by cavitation".

74-80, avril/mai.

70 Hydrodynamics - Concepts and Experiments

Process., 74 (2002) 177-187.

Fluid Engineering Centre, p.H1.1-H1.18.

Scientific Research (40); 1983. 377-391.

47-52.

ka, 1-8.

51-62.

Houille Blanche, Paris, année 47, n.316, juillet/aout.

