e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

available at www.sciencedirect.com

journal homepage: www.elsevier.com/locate/ecoleng

Modelling pollutant removal in a pilot-scale two-stage
subsurface ﬂow constructed wetlands
Attilio Toscano a,∗ , Günter Langergraber b , Simona Consoli a , Giuseppe L. Cirelli aa Department of Agricultural Engineering – Hydraulics division, University of Catania, Via S. Soﬁa 100, 95123 Catania, Italyb Institute of Sanitary Engineering and Water Pollution Control, BOKU – University of Natural Resources and Applied Life Sciences,
Vienna, Muthgasse 18, A-1190 Vienna, Austriaa r t i c l e

i n f oa b s t r a c t

Article history:

The hydraulic behaviour and efﬂuent pollutant concentrations in a pilot-scale two-stage

Received 8 December 2007

subsurface ﬂow constructed wetland for treatment of municipal wastewater have been sim-

Received in revised form 2 July 2008

ulated. The experimental pilot plant is located in San Michele di Ganzaria (Eastern Sicily)

Accepted 8 July 2008

and consists of four lines of two-stage subsurface ﬂow constructed wetlands for secondaryor tertiary treatment of municipal wastewater. The ﬁrst stage, for each line, consists of ahorizontal ﬂow bed, while in the second stage a vertical ﬂow bed operates for two lines

Keywords:

and a horizontal ﬂow bed for the other two. Phragmites sp. was used as vegetation in two

CW2D

lines while the other two lines are without plants. The HYDRUS-2D software was applied

Horizontal ﬂow

to describe ﬂow and single-solute transport, while the multi-component reactive transport

Numerical modelling

module CW2D was used to model the transformation and elimination processes of organicTwo-stage subsurface ﬂow

matter, nitrogen and phosphorus. Tracer studies and chemical wastewater analyses wereconstructed wetlands

carried out to calibrate and validate the transport model. In general, the simulation results

Vertical ﬂow

obtained show a good match with the measured data for water ﬂow, tracer experiments andpollutant removal processes.
© 2008 Elsevier B.V. All rights reserved.

1.

Introduction

Constructed wetlands (CWs) are by now recognised as effective and very suitable systems for natural (or extensive)
treatment of urban wastewater from small to medium communities (Kadlec et al., 2000) especially when the treatedwastewater shall be reused for crops irrigation (Barbagallo etal., 2001, 2002). Constructed wetlands are indeed characterised by low operation and maintenance costs and theiroperation requires only unskilled manpower. Furthermore,
the larger areas needed with respect to conventional wastewater treatment plants (WWTPs) usually do not represent a problem in regions, such as Southern Italy, with considerable landavailability in a rural context (Barbagallo et al., 2001). Although∗

Corresponding author. Tel.: +39 095 7147544; fax: +39 095 7147600.
E-mail address: attilio.toscano@unict.it (A. Toscano).
0925-8574/$ – see front matter © 2008 Elsevier B.V. All rights reserved.
doi:10.1016/j.ecoleng.2008.07.011

CWs are widely used and studied, they are often described as“black boxes” where the interactions between soil, vegetation,water and micro-organisms are not well known.
Usually CWs are designed and modelled using a simple
ﬁrst-order plug ﬂow reaction (PFR) model (USEPA, 1993; Kadlec
and Knight, 1996; Kadlec et al., 2000). As reported by several
authors (Shilton and Prasad, 1996; Bhattarai and Grifﬁn, 1999;
Kadlec, 2000), this approach ignores the considerable effectsof longitudinal dispersion, short-circuiting, dead zones androots resistance in the reaction volume of the bed. Instead of
the nominal hydraulic residence time (HRT) the actual ﬂow
conditions in the wetland can be better deﬁned by using theresidence time distribution (RTD) obtained by experimental
tracer tests (Werner and Kadlec, 2000; García et al., 2004).

282

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Fig. 1 – Layout of pilot plant.

According to Levenspiel (1972) and Kadlec and Knight(1996), the inﬂuence of the measured RTD can be consideredusing non-ideal ﬂow models such as “plug-ﬂow with dispersion” (PFD) or “tanks in series” (TIS). But these ﬁrst-order decay
models seem unable to accurately describe the ﬂow and theremoval processes occurring in the wetland, mainly due tothe strong interdependence between hydraulics and kinetics(Kadlec, 2000).
During the last few years, several numerical models withdifferent complexities have been developed as tools to get abetter understanding of the processes in CWs (Langergraber,2008). Only a few of these models are able to describe reactive
transport for variably saturated conditions that can be used formodelling vertical ﬂow as well as horizontal ﬂow in subsurface
constructed wetlands, and they have the highest complexity(Langergraber, 2001; McGechan et al., 2005; Freire et al., 2006;
Ojeda et al., 2006; Wanko et al., 2006).
In this paper the methodology proposed by Langergraber(2001) was used for simulating water ﬂow, tracer experiments
and removal processes in subsurface ﬂow constructed wetlands. The methodology was here applied to a pilot-scale
two-stage subsurface ﬂow constructed wetland, with horizontal and vertical ﬂow beds in series, located in Sicilyand founded by the Italian AQUATEC Project (Lopez et al.,2006). In particular, this method consists of the application of
the HYDRUS-2D/CW2D software (Langergraber and Šimůnek,2005): HYDRUS-2D (Šimůnek et al., 1999) is able to describe
ﬂow and single-solute transport in variably saturated porous
media, and the multi-component reactive transport module
CW2D (Langergraber, 2001, 2003) to model the transformation and elimination processes of organic matter, nitrogen andphosphorus in subsurface ﬂow constructed wetlands.

2.

Materials and methods

2.1.

Pilot plant

The experimental pilot plant used for this study is located in
San Michele di Ganzaria (Eastern Sicily), a rural community ofFig. 2 – Scheme of lines 1 and 2 of the pilot plant.

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

about 5000 inhabitants, and has been in operation since June
2004. Four parallel lines consist of two subsurface ﬂow beds inseries. The plant is designed to treat both efﬂuents from theprimary or secondary settler of the conventional wastewater
treatment plant (tricking ﬁlter).
Fig. 1 shows the layout of the pilot plant. Wastewater fromthe efﬂuent of the primary or secondary clariﬁer ﬁrst passes
a tank with a volume of 500 L and then is diverted to the fourlines. For each line the ﬁrst stage consists of a horizontal ﬂow(HF) bed, i.e. HF1-4, while in the second stage a vertical ﬂow(VF) bed operates for two lines (VF5-6) and a HF bed for theother two (HF7-8). Phragmites sp. with a density of 8 rhizomes
per m2 was used as vegetation in two lines (1 and 3) while theother two (2 and 4) are without vegetation.
The detailed scheme (plan, longitudinal and cross section)of lines 1 and 2 (horizontal + vertical systems) including theinlet and outlet structures is given in Fig. 2. Each bed has a rectangular shape with a surface area of 4.5 m2 (1.5 m × 3.0 m) andwas built in concrete and lined with an impermeable membrane. The horizontal ﬂow beds have a slope at the bottom of
about 1% and were ﬁlled, to an average depth of 0.6 m, with volcanic gravel having a size of about 10–15 mm. The water depth
in HF systems is about 0.55 m. The vertical ﬂow systems were
ﬁlled with two different medium layers: volcanic sand as main
layer for a depth of 0.5 m (about 0.06–4 mm; d10 = 0.13 mm;
d60 = 2.0 mm) at the top of the ﬁlter bed and coarse volcanic
gravel as drainage layer (about 16–32 mm, 25 cm depth).
Between the HF and VF beds in line 1 and 2, a tank collects
the efﬂuent of the HF bed for intermittent ﬂush loading of the
VF bed that is controlled by a water level sensor. The loadings
occur about every 2 to 2.5 h, a single loading amount 90 L. The
ﬁnal efﬂuent is used for irrigation of a green area close to thepilot plant.
The average inﬂuent ﬂow rate for each line was about
40 ± 5 L/h for the four HF beds (HF1-4) of the ﬁrst stage. At
the outlet of the experimental plant a probe and a data logger have been installed to record the electrical conductivity
and the temperature of efﬂuents from each line; a refrigerant
automatic sampler was also set-up at the inlet of the plant totake inﬂuent samples.
The research activity reported in this paper is mainly
focused on lines 1 and 2 (vegetated and not vegetated) functioning as tertiary treatment (ﬁrst experiment conducted) andsecondary treatment (second experiment). In particular, theﬁrst experiment was carried out just on line 2 using the pilot
plant as tertiary treatment, while the second experiment wason line 1 and 2 functioning as secondary treatment.

The ﬂow equation incorporates a sink term to account forwater uptake by plant roots. The solute transport equations
consider convective–dispersive transport in the liquid phase,
diffusion in the gaseous phase, as well as non-linear nonequilibrium reactions between the solid and liquid phases(Šimůnek et al., 1999).
The CW2D module considers 12 components and 9 processes. The components include dissolved oxygen, organic
matter (three fractions of different degradability, i.e. readily and slowly biodegradable, and inert), ammonium, nitrite,
nitrate, and nitrogen gas, inorganic phosphorus, and heterotrophic and two species of autotrophic micro-organisms.
Organic nitrogen and organic phosphorus are modelled as
nutrient contents of the organic matter, i.e. they are calculated as a percentage of COD. The biochemical elimination andtransformation processes are based on Monod-type expressions used to describe the process rates. All process rates anddiffusion coefﬁcients are temperature dependent.
The processes considered are hydrolysis, mineralization oforganic matter, nitriﬁcation (modelled as a two-step process),
denitriﬁcation, and a lysis process (as the sum of all decay
and loss processes) for the micro-organisms. CW2D assumes
a constant concentration of micro-organisms (and other compounds) in each ﬁnite element. The thickness of the bioﬁlm is
not considered.
The mathematical structure of CW2D is based on themathematical structure of the Activated Sludge Models, ASMs
(Henze et al., 2000). The total number of parameters in CW2D
is 46 (26 kinetic parameters, 6 stoichiometric parameters, 8
composition parameters describing the nitrogen and phosphorus content of the COD fractions, 2 parameters describing
the oxygen transfer between gaseous and aqueous phase,and 4 parameters describing the temperature dependency of
the kinetic parameters). Additionally, the material parametersand diffusion and adsorption parameters as used in HYDRUS2D have to be determined for each component. A detailed
description of the processes and parameters of CW2D is given
by Langergraber and Šimůnek (2005).
Langergraber (2005) investigated the plant uptake models
provided by HYDRUS that describe nutrient uptake coupledwith water uptake, and concluded that it was possible tosimulate plant uptake in high loaded systems, e.g. systems
treating mechanically pre-treated municipal wastewater. For
low-strength wastewater, the simulation results indicate that
potential nutrient uptake is overestimated by using thesemodels.

2.3.
2.2.

283

Measurements and calibration procedure

Simulation tool

The multi-component reactive transport module CW2D
(Langergraber, 2001) was developed to describe the biochemical transformation and degradation processes in subsurface
ﬂow constructed wetlands. CW2D was incorporated into the
HYDRUS-2D variably saturated water ﬂow and solute transport program (Langergraber and Šimůnek, 2005; Šimůnek etal., 1999).
The HYDRUS-2D program numerically solves the Richard’s
equation for saturated/unsaturated water ﬂow and theconvection–dispersion equation for heat and solute transport.

In order to calibrate and apply the described simulation tool,
several measurements and analyses have been conductedduring the investigation periods on line 1 and 2 (HF1 + VF5 and
HF2 + VF6, respectively) of the pilot plant: soil characteristics,
hydraulic properties of the media, ﬂow rate measurements,
tracer tests, physical and chemical analyses of wastewater.
To describe the hydraulic behaviour the following parameters
have been measured in the laboratory: porosity, bulk density
and saturated hydraulic conductivity.
For both the experiments, secondary and tertiary treatment, several pulse-input tracer tests were also conducted284

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

of a procedure proposed for the evaluation of the quality ofa model ﬁt. Further steps that have been applied and thatshowed a general good ﬁt are descriptive statistics and moresophisticated parameters such as the modiﬁed Coefﬁcient ofEfﬁciency.

3.

Results and discussion

3.1.

Numerical schematizations

The 2D ﬁnite element meshes used for the simulations areshown in Fig. 3. In particular, the horizontal ﬂow systems were
schematized by means of a longitudinal section (3.0 m lengthand 0.6 m depth) with 16 rows, 25 columns, 400 nodes and 720
ﬁnite elements (Fig. 3a), while the vertical ﬂow beds througha cross section (1.5 m width and 0.75 m depth) with 25 rows,
16 columns, 400 nodes and 720 ﬁnite elements (Fig. 3b). Atmospheric boundary conditions have been used for the inﬂuent
boundaries, constant pressure head boundary conditions forthe efﬂuent boundaries. The inﬂuent distribution in the vertical beds was restricted to two strips of 30 cm width each(Fig. 3b) to represent the inlet pipes of the vertical beds (seeFig. 2).

Fig. 3 – 2D mesh of HF1-2 (a) and VF5-6 (b).

using a potassium chloride (KCl) solution as tracer andmeasuring every 15 min the increases of electrical conductivity (EC) at the outlet of HF1-2 (by a portable conductivity
meter and data logger) and VF5-6 (by a ﬁxed probe), respectively.
Finally, HF1-2 inﬂuent and efﬂuent and VF5-6 efﬂuent were
sampled starting at the pilot plant start-up and the following
chemical-physical parameters were measured in laboratory
according to Standard methods (APHA, 1998): TSS (at 180 ◦ C),
BOD5 , COD, phosphorus (PO4 -P and TP), nitrogen (NH4 -N,
NO2 -N, NO3 -N and TN), EC, pH, dissolved oxygen (DO) andtemperature. To describe the dynamics of the system a few
intensive two-day measurement campaigns were planned andcarried out.
First, a calibration for tertiary treatment was carried out.
This could be only done for steady state conditions, as no
dynamics of the inﬂuent concentration occurred. As a second
step the calibration was performed for secondary treatment.
Again no dynamics of the inﬂuent concentrations occurred. To
test the dynamic behaviour of the model a pulse-input experiment for the NH4 -N inﬂuent concentration was performedby using an ammonium chloride (NH4 Cl) solution as reactive
tracer.
For estimating the match of the simulations to measured
data most of the times only a visual evaluation has been
carried out. According to Ahnert et al. (2007), visual evaluation (time series, scatterplot, plot of residuals) is the ﬁrst step

3.2.

Flow simulations

Table 1 shows the soil hydraulic parameters used for the ﬂow
simulations. The soil hydraulic parameters for the main layerof vertical beds in Table 1 were obtained using inverse simulations (Šimůnek et al., 1999) that use efﬂuent ﬂow ratemeasurements as input data.
The simulation results for water ﬂow obtained using theparameters given in Table 1 are shown in Fig. 4a using VF5 as
an example. Measured and simulated efﬂuent ﬂow rates of thevertical ﬂow bed VF5 are compared for two successive loadings. The simulated data show a good match with measured
data.
In Fig. 4b the overall ﬂow pattern of line 1 during one of theperformed tracer experiments is shown. The times of the loadings of VF5 have been recorded by a data logger. The inﬂuent
ﬂow rate for HF1 was calculated by the ﬁxed amount of a single
loading for VF5 (90 L) and the recorded intervals between two
successive loadings. Because of settling of suspended matterin the distribution pipe between two maintenance operations
the inﬂuent ﬂow rate for HF1 varied.

Table 1 – Soil hydraulic parameters for the van Genuchten modelParameter
Residual water content
r (dmw 3 /dms 3 )Saturated water content
s (dmw 3 /dms 3 )Saturated hydraulicconductivity Ks(m/h)

Empirical parameters

␣ (ms −1 )

n

HF1-2

0.050

0.175

500

2.78

4

VF5-6 – Main layer
VF5-6 – Drainage layer

0.045
0.050

0.30
0.19

0.36
500

2.78
2.78

3
4

285

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Table 2 – General transport model parameters for the tracer experiments
Parameter
Diffusion coefﬁcient(ms 2 /h)Longitudinal dispersion
coefﬁcient (m)Transversal dispersion
coefﬁcient (m)HF1-2

0.00050.45

0.1

VF5-6 – Main layer
VF5-6 – Drainage layer

0.00050.0005

0.15
0.15

0.5
0.1

Fig. 4 – Measured and simulated efﬂuent ﬂow rate of VF5 (a) and ﬂow pattern of line 1 (HF1 + VF5) (b).

3.3.
Single-solute transport simulations (tracer
experiments)
As mentioned before, several tracer tests were performed inlines 1 and 2. Table 2 reports the model parameters used for
simulating single-solute transport and Fig. 5 shows, as an
example, the results of some tracer experiments.
In particular, breakthrough curves of the measured andsimulated electrical conductivity in the vertical ﬂow beds VF5and VF6 are shown in Fig. 5a and Fig. 5b, respectively; Fig. 5cshows EC simulated and measured values for the horizontal
ﬂow bed at ﬁrst stage of line 2 (HF2), while the same data forthe overall two-stage constructed wetland (line 2) are reported
in Fig. 5d. The simulated breakthrough curves match the measured data well.
It is important to highlight that, to obtain the above
reported good results, the pore volumes of HF and VF main
layer used for the single-solute transport simulations had toFig. 5 – Tracer experiments for VF5 (a), VF6 (b), HF2 (c) and line 2 (HF2 + VF6) (d).

286e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Table 3 – Mean inﬂuent and efﬂuent concentrations measured and simulated during the experiments (standard
deviations in brackets)
Parameter(Unit)
Tertiary treatment

COD(mg/L)
Inﬂuent

NH4-N(mg/L)

NO2-N(mg/L)

NO3-N(mg/L)

Measured

22.9 (0.1)

8.6 (1.4)

0.83 (0.11)

8.4 (1.6)

Measured10.3 (1.9)

7.7 (0.3)

1.20 (0.09)

5.4 (1.3)

Simulated

10.0

8.96

0.03

7.0

Measured

8.4 (0.8)

4.9 (1.1)

1.00 (0.14)

20.0 (2.5)

Simulated

8.0

1.45

0.17

13.3

Measured

63.9 (3.3)

22.1 (1.9)

0.18 (0.38)

0.34 (0.10)

Measured

47.5 (2.8)

16.6 (1.2)

0.06 (0.05)

0.18 (0.17)

Simulated

44.4

17.2

0.02

2.82

Measured23.0 (2.5)

2.0 (0.1)

0.10 (0.07)

16.2 (1.1)

Simulated

24.6

0.54

0.01

26.9

Measured48.2 (4.1)

20.4 (0.6)

0.02 (0.01)

0.12 (0.11)

Simulated

50.6

21.5

6E-05

1E-04

Measured23.2 (2.3)

2.0 (0.1)

0.19 (0.08)

16.4 (1.3)

Simulated

24.0

0.96

0.01

25.0

Efﬂuent HF2 (unplanted)

Efﬂuent VF6 (unplanted)

Secondary treatment

Inﬂuent

Efﬂuent HF1 (planted)

Efﬂuent VF5 (planted)

Efﬂuent HF2 (unplanted)

Efﬂuent VF6 (unplanted)

be lower than measured pore volumes. In particular, for simulating the tracer experiments in the horizontal ﬂow systems(HF1-2) a signiﬁcant reduction of the porosity used for simulation with respect to that measured in laboratory (0.175 insteadof about 0.40) was needed. Obviously not all of the nominal
available reaction volume of the beds contributes to solutetransport and therefore to wastewater treatment. This could
be mainly due to phenomena like dead zones, preferential
ﬂow, hydraulic shortcuts and stagnant pore water that often
occur in subsurface wetlands but can not be easily measured
and described. This also implies that the hydraulic behaviourof the studied HF systems deviates considerably from ideal
plug ﬂow conditions as moreover reported by several authors(Kadlec and Knight, 1996; King et al., 1997; Kura et al.,
1997; Bhattarai and Grifﬁn, 1999; Barbagallo et al., 2003).
On the other hand, it is not possible that the horizontal
wetlands have actually lost more than the half of their porosity in the short time that they have been in operation (about
1 year between the plant start-up and the ﬁrst experimentconducted).

3.4.

Multiple-solute reactive transport simulations

The parameter set as proposed by Langergraber (2001) wasused for the reactive transport simulations. This parameter set
was developed and tested for the simulation of indoor pilotscale subsurface vertical ﬂow CWs with intermittent loading
treating mechanically pre-treated wastewater (Langergraber,
2003).
Table 3 shows the mean measured and simulated inﬂuent and efﬂuent concentrations of lines 1 and 2 during thedifferent 2-day sampling periods for tertiary and secondary
treatment experiments, respectively. The water temperatureduring the two experiments was 6.5 ◦ C (tertiary) and 24.1 ◦ C(secondary treatment). Low inﬂuent concentrations have been
detected probably because the raw wastewater was diluted
with inﬁltrating water coming from a close aqueduct during
the measurement campaign.
To estimate the organic matter fractions for the model
input (CR = ready biodegradable COD, CS = slowly biodegradable COD and CI = inert COD), the following relations havebeen used for inﬂuent fractionation: CI = 0.85 × CODefﬂuent ;CR = 2 × CS; COD = CR + CS + CI.
Quasi-steady-state conditions have been obtained after
about three months of simulation time for each experiment. A
constant inﬂuent ﬂow rate was used for the HF beds except for
the last 2 weeks prior the sampling campaign. For this period
the recorded times for loading of vertical systems have beenused to calculate the inﬂow of HF1-2. This results in a variable
ﬂow rate while the concentrations remained constant.
For tertiary treatments, simulation results showed a good
match between measured and simulated data for COD concentrations for both HF2 and VF6 efﬂuents. Nitrogen simulations
seem to be less accurate. Ammonia removal was slightly
underestimated for HF2 and signiﬁcantly overestimated for
VF6, while the simulation of nitrate has a good ﬁt for HF2 andis underestimated for VF6. This can probably be explainedby the fact that the standard kinetic parameters used havebeen previously tested for secondary treatment only, while
the pilot-scale plant has operated as tertiary treatment.
Additionally, the Activated Sludge Models are deﬁned forwater temperature higher than 10 ◦ C (Henze et al., 2000).
However, the main measurements used for simulation
were instead conducted during the winter period when
the wastewater temperature was outside this temperature
range.

287

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Table 4 – Simulated evapotranspiration, O2 release and nitrogen uptake
Parameter
Evapotranspiration

HF1
VF5

(L/d)

(L/m2 /d)

33.1
33.3

7.4
7.4

O2 release

NH4 uptake(g/d)(g/m2 /d)

22.5
22.7

5.01
5.04

Fig. 6 – Measured inﬂuent concentrations and measured
and simulated HF2 and VF6 efﬂuent concentrations of COD(without vegetation, secondary treatment).

For secondary treatment the simulations were carried out
for line 1 and 2, i.e. with and without vegetation, respectively.
The effects of plants were considered by taking in accountoxygen release and nutrients uptake. The root zone in the HF1
bed was set to be 30 cm deep, i.e. half of the depth of the horizontal system (Rousseau, 2005). O2 release was assumed tobe 5 g O2 /m2 /d (Vymazal et al., 1998). The evapotranspiration
rate was calculated to be 7.4 mm/d for the investigation period(July) using a modiﬁed Penman–Monteith equation (ASCEEWRI, 2005) with the meteorological data recorded througha weather automatic station installed in the pilot plant area.
For the simulations it was assumed that ammonia and nitrate
are taken up by plants with the same preference. In HYDRUS2D uptake and release of solutes by plants is coupled to water
uptake (Šimůnek et al., 1999). O2 release was therefore calibrated to reach the target value (Table 4). The simulation
predicts that nitrogen is mainly taken up as ammonia in HF1
bed and as nitrate in the VF5, respectively.
Figs. 6 and 7 show measured and simulated COD and NH4 N concentrations for line 2 without vegetation (HF2 + VF6) of
the pilot plant as secondary treatment. The inﬂuent concentrations showed almost no variation. The graphs highlight avery good match between measured and simulated values of
COD concentrations (Fig. 6). Also simulation results for NH4 N seem to be quite accurate, almost no NH4 -N removal for
HF2 and slightly overestimated NH4 -N removal for VF6 (Fig. 7).
The match of simulated ammonia to measured data was better as compared to tertiary treatment. This was expected
since the investigations concerning secondary treatment havebeen carried out in summer with water temperature around24 ◦ C.

NO3 uptake(g/d)(g/m2 /d)(g/d)(g/m2 /d)

0.68
0.04

0.15
0.01

0.05
0.78

0.01
0.17

Fig. 7 – Measured inﬂuent concentrations and measured
and simulated HF2 and VF6 efﬂuent concentrations of
NH4 -N (without vegetation, secondary treatment).

Figs. 8 and 9 show measured inﬂuent concentrations andmeasured and simulated efﬂuent concentrations for HF1 andVF5 (line 1, with vegetation) for COD and NH4 -N, respectively. In both ﬁgures two simulation results are shown: withand without considering the effects of plants in the simulation. Simulated COD efﬂuent concentrations for HF1 areslightly lower when considering plants in the simulation(Fig. 8). Ammonia efﬂuent concentrations for HF1 appear signiﬁcantly overestimated (no removal at all) in the simulation
if the effects of the plants are not considered. The measured
data can only be matched when considering oxygen releaseand nutrient uptake (Fig. 9). Simulated NH4 -N efﬂuent concentrations for VF5 are slightly underestimated for both thesimulations.

Fig. 8 – Measured inﬂuent concentrations and measured
and simulated HF1 and VF5 efﬂuent concentrations of COD
without and with plants (secondary treatment).

288

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Fig. 9 – Measured inﬂuent concentrations and measured
and simulated HF1 and VF5 efﬂuent concentrations of
NH4 -N without and with plants (secondary treatment).

Comparing the simulation results obtained for line 1 andline 2, an inﬂuence of vegetation in the ammonia removal
is evident for HF beds. A moderate NH4 -N removal wasmeasured in HF1 (vegetated) and could be simulated considering release of oxygen and nitrogen uptake by plants. No
ammonia removal was observed and simulated in HF2 (notvegetated) where anaerobic conditions fully occur due to theabsence of plants (measured values of dissolved oxygen in
the bed always < 1.0 mg/L). Most of the ammonia nitrogen wasremoved in the second VF stage of both the studied lines bynitriﬁcation.
Conversely, COD removal was only slightly inﬂuenced by
the presence of vegetation (such as reported by e.g. Tanner,2001) and occurred mostly in the ﬁrst horizontal stage during the tertiary treatment experiment, while it was detected
mainly in the second stage for secondary treatment. This wasprobably due to the very low inﬂuent COD values in the tertiary
treatment that had already been removed to the “backgroundconcentration” in the ﬁrst stage.
Since no signiﬁcant variation occurred for in the inﬂuent
pollutant concentrations, a further experiment was conductedsimulating inﬂuent ammonia variations and checking the

Fig. 10 – Measured inﬂuent concentrations and measured
and simulated efﬂuent concentrations of NH4 -N after the
addition of NH4 Cl solution in line 1 (secondary treatment).

Fig. 11 – Measured inﬂuent concentrations and measured
and simulated efﬂuent concentrations of NH4 -N after the
addition of NH4 Cl solution in line 2 (secondary treatment).

results of the model. This was done by adding ammonium
chloride (NH4 Cl) as pulse tracer in the inﬂuent of lines 1
and 2 functioning as secondary treatment. During this period
the hydraulic loading of the HF1 was 900 L/d, i.e. 200 mm/d,
about twice as during the tracer experiments described
before.
Measured and simulated NH4 -N concentrations arereported for line 1 (vegetated) in Fig. 10 and for line 2(unplanted) in Fig. 11. The previous values for the parameters describing plant uptake and oxygen release (Table 4) havebeen used for the simulations. CW2D seems to be able to predict the dynamics in the efﬂuent of the HF bed for both lines. In
particular, it was possible to match measured values of HF1-2
efﬂuent also during the peak after the addition of solution. As
observed before, a small ammonia removal occurred only in
HF1 (vegetated) while the most of NH4 -N reduction occurredin the vertical beds.

4.

Conclusion

A simulation study of ﬂow, solute transport and removal
processes in pilot-scale two-stage subsurface constructedwetlands used for secondary/tertiary treatment of urban
wastewater resulted in good correlations with measured
indices of hydraulic behaviour and biochemical eliminationand transformation processes. These results conﬁrm thatthe hydraulic behaviour of the system must be correctly
described to obtain a good match between modelled andmeasured efﬂuent parameters. In particular, efﬂuent COD
concentrations were always close to measured values bothfor the secondary and tertiary treatment. Ammonia simulations pointed out unsatisfactory results during tertiary
treatment, while in the secondary treatment experiment NH4 N efﬂuents values were simulated in an adequate manner.
When considering the effects of oxygen release and nutrient
uptake by roots, the differences in the NH4 -N efﬂuent concentrations in planted/unplanted HF beds can be simulated.
The model seemed also able to simulate ammonia efﬂuent
dynamic.

e c o l o g i c a l e n g i n e e r i n g 3 5 ( 2 0 0 9 ) 281–289

Acknowledgements
This research was funded, within the AQUATEC Project, by
the Italian Ministry of Instruction, University and Research(MIUR) and by the European Union. The Authors thank Mr
Mirco Milani, Dr Vanessa Di Grande and Ms Alessia Marzo
for carrying out the experiments and for analysing thesamples.

references

Ahnert, M., Blumensaat, F., Langergraber, G., Alex, J., Woerner, D.,
Frehmann, T., Halft, N., Hobus, I., Plattes, M., Spering, V.,
Winkler, S., 2007. Goodness-of-ﬁt measures for numericalmodelling in urban water management—a review to support
practical applications. In: Proceedings of the 10th IWA
Specialised Conference on Design, Operation and Economicsof Large Wastewater Treatment Plants, 9–13 September 2007,
Vienna, Austria, pp. 69–72.
APHA, 1998. Standard Methods for the Examination of Water andWastewater, 20th ed. American Public Health Association,
Washington D.C., USA.
ASCE-EWRI, 2005. The ASCE Standardized reference
evapotranspiration equation. ASCE-EWRI Standardization ofreference evapotranspiration. Task committee Final Report,
available online at http://www.kimberly.uidaho.edu/
water/asceewri/.
Barbagallo, S., Cirelli, G.L., Indelicato, S., 2001. Wastewater reusein Italy. Water Sci. Technol. 43 (10), 43–50.
Barbagallo, S., Cirelli, G.L., Consoli, S., Toscano, A., Zimbone, S.M.,
2002. Performances of a H-SSF constructed wetland as tertiary
treatment for wastewater reuse: the case study of “S. Michele
di Ganzaria” in Sicily, Italy. In: Proceedings of IWA Regional
Symposium on Water Recycling in Mediterranean Region,
26–29 September, Heraklion, Greece,
pp. 379–385.
Barbagallo, S., Cirelli, G.L., Consoli, S., Toscano, A., Zimbone, S.M.,
2003. Performances and Hydraulics of a H-SSF constructedwetland for municipal wastewater reuse. In: Proceedings of
4th IWA International Symposium on Wastewater
Reclamation and Reuse, 12–14 November, Mexico City, Mexico.
Bhattarai, R.R., Grifﬁn, D.M., 1999. Results of tracer tests inrock-plant ﬁlters. J. Environ. Eng. 125 (2), 117–125.
Freire, F.G., Davies, L.C., Vacas, A., Pedro, I., Novais, J.M.,
Martins-Dias, S., 2006. Continuous and intermittent loading ofa vertical ﬂow constructed wetland for an azo dye treatment.
In: Proceedings of the 10th IWA Specialized Group Conference
on “Wetland Systems for Water Pollution Control”, vol. 3,
23–29 September, Lisbon, Portugal,
pp. 1501–1509.
García, J., Chiva, J., Aguirre, P., Alvarez, E., Sierra, J.P., Mujeriego,
R., 2004. Hydraulic behaviour of horizontal subsurface ﬂow
constructed wetlands with different aspect ratio and granularmedium size. Ecol. Eng. 23, 177–187.
Henze, M., Gujer, W., Mino, T., van Loosdrecht, M.C.M., 2000.
Activated sludge models ASM1, ASM2, ASM2D and ASM3. IWA
Scientiﬁc and Technical Report No. 9. IWA Publishing, London,UK.
Kadlec, R.H., 2000. The inadequacy of ﬁrst-order treatment
wetland models. Ecol. Eng. 15, 105–119.
Kadlec, R.H., Knight, R.L., 1996. Treatment wetlands. Lewis
Publishers, New York, USA.

289

Kadlec, R.H., Knight, R.L, Vymazal, J., Brix, H., Cooper, P., Haberl,R., 2000. Constructed Wetlands for Pollution Control. IWA
Scientiﬁc and Technical Report No. 8. IWA Publishing, London,
UK.
King, A.C., Mitchell, C.A., Howes, T., 1997. Hydraulic tracer studies
in a pilot scale subsurface ﬂow constructed wetland. Water
Sci. Technol. 35 (5), 189–196.
Kura, B., La Motta, E.J., Tittlebaum, M., Alawady, M., 1997.
Macroscopic BOD kinetic model for microbial rock plant ﬁlter
design. Adv. Environ. Res. 1 (1), 36–43.
Langergraber, G., 2001. Development of a simulation tool forsubsurface ﬂow constructed wetlands. Wiener Mitteilungen169. Vienna, 207 pp, ISBN 3-85234-060-8.
Langergraber, G., 2003. Simulation of subsurface ﬂow constructedwetlands—results and further research needs. Water Sci.
Technol. 48 (5), 157–166.
Langergraber, G., 2005. The role of plant uptake on the removal oforganic matter and nutrients in subsurface ﬂow constructedwetlands—a simulation study. Water Sci. Technol. 51 (9),213–223.
Langergraber, G., 2008. Modeling of processes in subsurface ﬂow
constructed wetlands—a review. Vadoze Zone J. 7 (2), 830–842.
Langergraber, G., Šimůnek, J., 2005. Modeling variably saturated
water ﬂow and multicomponent reactive transport inconstructed wetlands. Vadoze Zone J. 4 (4), 924–938.
Levenspiel, O., 1972. Chemical Reaction Engineering, second ed.
Wiley, New York, USA.
Lopez, A., Pollice, A., Lonigro, A., Masi, S., Palese, A.M., Cirelli,G.L., Toscano, A., Passino, R., 2006. Agricultural wastewaterreuse in southern Italy. Desalination 187, 323–334.
McGechan, M.B., Moir, S.E., Castle, K., Smit, I.P.J., 2005. Modelling
oxygen transport in a reedbed-constructed wetland
puriﬁcation system for dilute efﬂuents. Biosyst. Eng. 91 (2),
191–200.
Ojeda, E., Caldentey, J., García, J., 2006. 2D simulation model forevaluating biogeochemical pathways involved in organic
matter removal in horizontal subsurface ﬂow constructedwetlands. In: Proceedings of the 10th IWA Specialized Group
Conference on “Wetland Systems for Water Pollution Control”,
vol. 3, 23–29 September, Lisbon, Portugal, pp. 1405–1413.
Rousseau, D.P.L., 2005. Performance of constructed treatment
wetlands: model-based evaluation and impact of operationand maintenance. Ph.D. Thesis. Ghent University, Belgium.
Shilton, A.N., Prasad, J.N., 1996. Tracer studies of a gravel bed
wetland. Water Sci. Technol. 34 (3–4), 421–425.
Šimůnek, J., Senja, M., van Genuchten, M.Th., 1999. The
HYDRUS-2D software package for simulating thetwo-dimensional movement of water, heat, and multiple
solutes in variably-saturated media. Version 2.0 Manual. U.S.
Salinity Laboratory, USDA, Riverside, California, USA.
Tanner, C.C., 2001. Plants as ecosystem engineers insubsurface-ﬂow treatment wetlands. Water Sci. Technol. 44(11–12), 9–17.
USEPA, 1993. Subsurface Flow Constructed Wetlands for
Wastewater Treatment—A Technology Assessment.
EPA/832/R-93/008. U.S. Environmental Protection Agency,
Washington DC, USA.
Vymazal, J., Brix, H., Cooper, P.F., Green, M.B., Haberl, R. (Eds.),
1998. Constructed wetlands for wastewater treatment inEurope. Backhuys Publishers, Leiden, The Netherlands.
Wanko, A., Mose, R., Carrayrou, J., Sadowski, A.G., 2006.
Simulation of biodegradation in inﬁltration seepage—model
development and hydrodynamic calibration. Water Air Soil
Poll. 177, 19–43.
Werner, T.M., Kadlec, R.H., 2000. Wetland residence time
distribution modelling. Ecol. Eng. 15, 77–90.