Influence factors in the adjustment of parameters of the modified first-order kinetics equation used to model constructed wetland systems

The aim of the present work was to evaluate the relationship between the hydraulic retention time power parameter (n) of the modified first-order kinetic mathematical model (MFOKM) with the hydraulic dispersion number (d) and the organic matter degradability, obtained from the ratio of chemical and biochemical oxygen demands (COD/BOD) along the length of horizontal subsurface flow constructed wetlands (HSSF-CWs). For this purpose, unplanted HSSF-CWs prototypes of equal volumes and length/width ratios (L W-1) of 1.0, 4.0 and 7.3 were used. All HSSF-CWs were filled with gravel # 0 and received sanitary sewage (SW) during a period of 90 days, in which the BOD and COD, in the total and soluble forms, were monitored in the influent and effluents from the system. These units were subject to tracer tests, before and after the HSSF-CWs received the SW, to obtain d. The MFOKM was adjusted to the COD and BOD experimental data. In all cases, an inverse relation was observed between n and the parameters d or with the organic matter decay rate obtained in the initial part of the HSSF-CWs.


Introduction
Constructed wetlands systems (CWs) have been presented as a potential alternative for the treatment of wastewaters.Among the various types of constructed wetlands systems (CWs), those with horizontal subsurface flow (HSSF-CW) are the most common in the USA and Europe and studied in Brazil (Machado, Beretta, Fragoso, & Duarte, 2017).These systems differ from the other configurations by the fact that the wastewater flows through the support medium horizontally, from inlet to outlet, without exposure to the atmosphere.
The design of such systems is based on empiricism (Langergraber et al., 2009), notably in relation to the criteria of required area, loading or admissible per capita ratios, or the pollutant decay equations (Kadlec, & Wallace, 2008).
The First Order Kinetic Model (FOKM) is the most utilized model for describing the process (Kadlec & Wallace, 2008), however it has not shown good fit to the experimental data because the reactors do not presented an ideal hydraulic behaviour, due to the occurrence of flow anomalies such as dispersion, shortcircuits and the presence of dead zones (Dittrich & Klincsik, 2015), and due to the decay coefficient of organic material suspended in wastewater which varies with the hydraulic retention time (HRT) (Shepherd, Tchobanoglous, & Grismer, 2001).According to Langergraber et al. (2009), the simplified assumption of an ideal plug-flow in HSSF-CWs may introduce significant errors to the parameters, as well as the HRT, and consequently on the pollutant removal estimates.
It has recently been suggested that the interdependence between hydraulics and kinetics is very strong and influenced by several factors, and that even the latest generation of mathematical models has been considered inadequate to fully explain the observed behaviour and provide a reliable estimate of the parameters (Marsili-Libelli & Checchi, 2005, Matos, Matos, Costa, & von Sperling, 2018).Therefore, it is necessary to expand the knowledge regarding the hydraulic behaviour of HSSF-CWs in order to improve estimates of the kinetic removal rate and to obtain satisfactory predictions of pollutant efficiency in these systems (Garcia et al., 2004).For these reasons, some adaptations to the FOKM have been proposed (Laber, Habert, & Shrestha, 1999, Shepherd et al., 2001, Brasil, Matos, Silva, Cecon, & Soares, 2007, Matos et al., 2018).
According to Brasil et al. (2007) and Matos et al. (2018), variation in the decay coefficient of the organic material in the wastewater under treatment in HSSF-CWs can be associated to both the decrease in the easily degradable organic fraction, which occurs as the wastewater moves towards the outlet of the reactor, and to changes in the hydraulic behaviour inside the system, especially with regards to the decrease in the HRT.These authors proposed the insertion of a coefficient n, as a power parameter of the HRT variable, in the FOKM, to give Equation 1 which expresses the Modified First Order Kinetic Model (MFOKM) and provides good fit to the experimental data obtained by different authors and several different wastewaters (Brasil et al., 2007, Matos et al., 2018).

= .
. (1) The coefficient n, according to the authors, has the function of correcting the inadequacies arising from the reduction in the biodegradability of the organic material along the length of the HSSF-CW, as well as the dispersive effects caused by the flow, through an adjustment of the HRT.
According Matos et al. (2018), based on the suggestion of von Sperling by personal communication in 2016, the analysed MFOKM equation (Equation 1) developed by Brasil et al. (2007) proposed developments of the exponential term, in a similar way to what was undertaken by Shepherd et al. (2001), i.e. the insertion of an adjustment for the organic matter decay coefficient with the HRT, so that it would be possible to show the cited model in the FOKM format, as shown in Equation 2.
where k B = MFOKM decay coefficient (d -1 ); k f(HRT) = variable decay coefficient as a function of HRT, to be used in the corrected FOKM (d -1 ).
Thus, k f(HRT) is the decay coefficient, varying with the HRT of the system, which incorporates the decrease in biodegradability of the remaining organic matter in the wastewater.With the use of k f(HRT) , the HRT exponent can be maintained equal to 1, thus eliminating the variable n and maintaining the original FOKM equation.
It is known that the length-width ratio (L W -1 ) has a strong influence on the hydraulic behaviour of a reactor, and should therefore influence the adjustment of kinetic models developed to describe the removal of organic matter when used for wastewater treatment (von Sperling, 2007).However, few studies have addressed the influence of reactor configurations in this regard.Alcocer, Vallejos, and Champagne (2012) found a large influence of the L W -1 ratio on the hydraulic behaviour of HSSF-CWs prototypes, stating that an increase in this ratio provides behaviour more similar to plug-flow.On the other hand, Garcia et al. (2005) and Zhang et al. (2014) concluded that the HRT and CW depth are factors of greater influence on the efficiency of pollutant removal in these systems than the L W -1 .Ranieri, Gorgoglione, and Solimeno (2013) observed that unplanted HSSF-CWs maintain a plug-flow regime, however in similar systems cultivated with Phragmites and Typha, the plants provide a deviation from this flow regime.
Hydrodynamic characterization of the CWs was evaluated by performing tracer tests, and after analyzing the breakthrough curves (BTCs) of the effluents, the hydraulic retention time, porosity, variance and dispersion number parameters can be obtained, which may reveal the current system's operating condition, the presence of short circuits or dead zones, and allow for comparisons between different designs and methods of operating these systems (Seeger, Maier, Grathwohl, Kuschk, & Kaestner, 2013).
According to Garcia et al. (2005), in HSSF-CWs the d varies with the hydraulic loading rate, HRT, L W -1 ratio, cultivated macrophyte type, wastewater quality, system operating time, clogging degree and others.The flow dispersion number, represented by the parameter 'd', is dimensionless and is considered intermediate if greater than 0.025 and high if greater than 0.2, indicating a high degree of mixing (Kadlec & Wallace, 2008).
Therefore, the objective of this work was to investigate the association of the variable n with the hydrodynamic conditions of the reactor, mainly with the dispersion number d and also with the decrease in the degradability of the organic material along the length of three horizontal subsurface flow constructed wetland systems with different geometric configurations.

Material and methods
The experiment was carried out in a greenhouse located at the Department of Agricultural Engineering of the Federal University of Viçosa (UFV), Viçosa, State of Minas Gerais, Brazil, at the following geographic coordinates: latitude 20° 45' 14" S and longitude 42° 52' 53" W and average elevation of 650 m above sea level.
The influent used in the study was raw sewage (SW) drawn from the sewerage network of the Bosque Acamari condominium and its surroundings, whose main characteristics are show in Table 1.The SW was pumped to a storage reservoir, from which it was pumped into a 1000 L reservoir, installed inside the greenhouse, every 5 days.From this reservoir, it was distributed to the HSSF-CWs through three ProMinet ® CONCEPT peristaltic pumps, operated at an average flow rate of 0.03 m 3 d -1 .This flow rate provided a nominal HRT (HRT n ) in the CW units of approximately 3 days, a period considered sufficient to evaluate the organic material removal process in the evaluated systems.
The HSSF-CWs used in the evaluations were constructed on a pilot scale, consisting of fiberglass boxes, all with a total volume of 0.40 m 3 and a height of 0.40 m.The length and width dimensions varied among the HSSF-CWs, being 1.00 x 1.00 m in CW1, 2.00 x 0.50 m in CW2, and 2.70 x 0.37 m in CW3, which provided length-to-width (L W -1 ) ratios of 1.0, 4.0 and 7.3, respectively, as shown in Table 2.
The HSSF-CWs were filled with a 0.25 m layer of pea gravel and the water level was maintained at 0.20 m, therefore leaving a freeboard of 0.05 m.The support medium used to fill all CWs consisted of gravel # 0 (diameters D 10 = 4.47 mm, D 60 = 9.25 mm uniformity coefficient -CU D 60 /D 10 = 2.07 and initial void or pore volume of 0.420 m 3 m -3 ), reused from previously conducted experiments with HSSF-CWs.Before being added to the fiberglass boxes, the gravel was passed through a 2 mm mesh sieve to remove the fine fraction.The units were maintained unplanted so that there was no other factor influencing the hydraulic characteristics and degradation of the organic matter (Matos, von Sperling, Matos, & Passos, 2015;Matos et al., 2018).The HSSF-CWs operated in parallel, where the influent was distributed by a PVC tube (25 mm) which entered the units at half the liquid column depth (0.10 m from the bottom), as shown on Figure 1.At that position, a 'T' connection was installed in which a perforated tube was fitted which distributed the SW along the full width at the inlet of the HSSF-CWs.The effluent was collected from the bottom of the opposite end, also by a perforated tube arranged to cover the full width.On the outside, a vertical tube was coupled to for SW level control inside the HSSF-CWs.
A grid of collection points was installed along each HSSF-CW, consisting of 32 mm diameter PVC tubes, which had a perforated base to allow free entry of the SW into the sample collection tube.These tubes were arranged 2 by 2, at four equidistant positions along the direction of flow in the HSSF-CWs (Figure 1).
The experiment was conducted between June and September of 2015 when the SW application was initiated, with the tracer test performed between August and September.
The methodology used for the hydrodynamic tracer test was established by Headley and Kadlec (2007), and consisted of using lithium chloride salt (LiCl) as a tracer, with pulse injection for approximately 2 min.in each system.The tracer solution was prepared at the concentration of 166.7 g L -1 of LiCl, established according to the background concentration of the influent, the detection limit curve of the equipment and the desired maximum concentration peak of the tracer (Headley & Kadlec, 2007).As soon as the tracer solution was injected, samples were collected at pre-established intervals and stored in small plastic bottles for subsequent reading of the Li concentration in the samples, in a flame photometer (CELM-180 ® ).Collections were performed at the following sample frequency: on the 1st day of collection sampling every 1 hour, from the 2 nd to 6 th day sampling every 2 hours, from the 7 to 11 th day every 4 hours, from the 12 to 15 th days every 8 hours, and from the 16th to 21st day every 24 hours.The data obtained were considered until the LiCl concentration in the effluent became negligible.
From the concentration data and the permanence time of each sample the C pulse curve was constructed, from which the measured hydraulic retention time (HRT m ) and the variance (σ 2 ) were extracted, used to calculate the dispersion number (d) according to Equation 3, presented by Levenspiel (2011).
HRT mr 2 = 2d + 3d 2 (3) The boundary condition considered in the mathematical solution that generated Equation 3, used in the calculation of d, was that of an open-closed vessel, where the tracer solution was added to a tubing that leads directly to the interior of the support medium and the effluent was collected in the outlet tube, outside the CW.
It is important to emphasize that the support medium used had already been used in other SW experiments, and therefore already had an attached biofilm, which accelerated the biomass stabilization process.Sampling and analysis campaigns during the system monitoring period occurred in the months of July and August, where six samples were always collected in the morning, with a frequency of 10 days.During these campaigns, samples of the influent, effluent and intermediate positions of the HSSF-CWs were collected.In the case of the influent, a single sample was taken in each campaign, since it consisted of wastewater from the same source.Because each intermediate sampling position was formed by two parallel tubes, the representative sample of each collection position was composed of 300 mL from each.The effluent samples were collected at the exit of each HSSF-CW.
The collected samples were packaged in an icebox and then taken to the Water Quality Lab of the Agricultural Engineering Department, UFV, for analysis.Characteristics including temperature and flow were measured in loco, where the temperature was measured by means of a mercury thermometer and flow by the direct volumetric method.
In the Water Quality Lab, in accordance with the recommendations of the Standard Methods for Water and Wastewater Examination (American Public Health Association/American Water Works Association/Water Environment Federation [Apha/Awwa/WEF], 2012), analyses were conducted of total and soluble COD by the open reflux chemical oxidation method, and total and soluble BOD 5 by the Winkler method.The total and soluble COD and BOD variables were analyzed considering water losses by evaporation, assuming that these happened homogeneously in all HSSF-CWs, with the concentration corrections made proportionally to the distance of each sampled position, as shown on Equation 4.
where: Q 0 = influent flow to the system (m 3 d -1 ); Q e = effluent flow from the system (m 3 d -1 ), c = distance from inlet to any longitudinal position inside the units (m), L = total length of the unit (m), Cc = corrected concentration for water losses (mg L -1 ) and; Cm = measured concentration (mg L -1 ).
For the evaluation of intensity of the correlation between variables, the Pearson coefficient was used.The influence of the L W -1 ratio on the adjusted parameters was evaluated using the Kruskal-Wallis.

Results and discussion
The d values obtained in this work, for the open-closed contour condition, were 0.13, 0.10 and 0.09, respectively, for CW1, CW2 and CW3 and were within the moderate dispersion range (0.025 < d < 0.20), as established by Kadlec and Wallace (2008).
Influent and effluent temperatures of the HSSF-CWs were in the range of 19 to 26 ºC during the experimental period, and it was verified that the geometric configuration of these reactors did not have a sensible effect on this variable.
The effluent flow from the HSSF-CWs was measured on the days of sample collection, totalling 6 measurements with the following means (standard deviations): 0.0281 (0.0018) m 3 d -1 in CW1, 0.0288 (0.0024) m 3 d -1 in CW2 and 0.0297 (0.0032) m 3 d -1 in CW3.Water losses by evaporation, measured by the difference between the influent and effluent flows, were 7.1 in CW1, 4.8 in CW2 and 1.8% in CW3.The lowest water losses were observed in the unit with the elongated reactor (greater L W -1 ratio), despite presenting the same exposed surface area.Although they are relatively small, these evaporative water losses were used, as presented in Equation 5, to correct the total and soluble COD and BOD along the reactors.
It is verified that all HSSF-CWs showed a similar behaviour in terms of organic matter decay along the length of the reactor, and it is not possible to identify an influence of the HSSF-CWs geometric configuration.This somewhat goes against what is established in terms of system performance under its geometric configuration, which conditions the hydraulic drainage model.Based on known theories on the subject, reactors with greater L W -1 ratios should approach a plug-flow regime, thus leading to a greater decay of the organic matter, represented by a first-order reaction.The lower the L W -1 ratios, the closer the hydraulic regime is to complete mixing, leading to a lower decay of the pollutants.However, the behaviour observed in this work is corroborated by Garcia et al. (2005), who stated that because of organic matter removal in HSSF-CWs occurring mainly through physical mechanisms (sedimentation and filtration), the removal efficiency of this pollutant is very similar between systems with different L W -1 , when maintaining the same water depth in the systems.
The MFOKM was adjusted to the data of each analysed variable and presented in the respective HSSF-CW curves, and the parameters are shown in Table 3.
Comparing the fitting of the MFOKM to the experimental data by means of the Coefficient of Determination (R 2 ), a good fit is verified, providing values always superior to 0.98.In terms of the k B and n parameters adjusted in the MFOKM, the values were maintained between 0.52 and 0.74 d -1 , close to those obtained by Chagas, Matos, Cecon, Monaco, and França (2011) which were 0.45 and 0.69 d -1 when treating effluents from a septic tank in HSSF-CWs under subsurface application rates (TAS) of 44 and 98 kg ha -1 d -1 of COD.
Although there was no statistical difference between data according to the Kruskal Wallis test, a tendency for k B to decrease and n to increase with increase in the L W -1 ratio was observed.The same was verified by Matos et al. (2018) when adjusting the degradation coefficient of the Residual Concentration Model (RCM) equation, observing a decrease of the parameter value with the increase of the L W -1 ratio.Because k B and n individually express the influence on decay of organic matter in the residual water under treatment, as suggested by Brasil et al. (2007) and can be seen in Equation 1, the joint influence of both parameters is what should be evaluated, since it expresses the parameter k f(HRT) of Equation 2.
Regarding the values of the n parameter adjusted in the MFOKM to the COD and BOD data, the values remained at around 0.42 to 0.83.Brasil et al. (2007), when studying HSSF-CWs cultivated with a local reed and L W -1 equal to 24, in the treatment of effluents from a sceptic tank, values of n were found between 0.20 and 0.25.Chagas et al. (2011), in HSSF-CWs with the same L W -1 ratio, but cultivated with yellow lily and submitted to different organic loads, obtained values of 0.11 to 0.46.
Analysing the values adjusted for the parameter n, its relatively small variability can be verified (values ranging from 0.43 and 0.70) referring to the curves adjusted for BOD decay, both raw and soluble, and with variability slightly higher than adjustment of the COD curves.In both cases, the results obtained indicated that the CW geometry did not permit for a significant change or did not show identifiable association with the n parameter value.Note: CW1 square reactor (L W -1 = 1.0);CW2 intermediate rectangular reactor (L W -1 = 4.0); CW3 elongated reactor (L W -1 = 7.3).
To determine if the degradability of organic matter along the HSSF-CWs may be a factor of influence on the n value, the COD/BOD ratio was evaluated as a representative variable of such degradability.The higher this relationship, the higher the recalcitrance of the organic materials present in residual water under treatment.Figure 2 presents boxplots representing this relationship along the length of the HSSF-CWs.
Although elevated standard deviations were observed in all sample points, the behaviour of this ratio shows a reasonably stability in all graphics, without presenting decay of growth tendency along the length of HSSF-CWs, which is an indication that the decay of these two variables occurs at similar magnitudes.Another evaluation performed indicated that there was no sensitive variation in biodegradability of SW under treatment along the length of the HSSF-CWs.This result is different from what was expected, since the more labile fraction of organic matter (expressed by BOD) should be removed at higher rates, notably at the beginning of the HSSF-CWs (Matos, Freitas, & Borges, 2011).In literature, a reduction in the COD/BOD ratio can also be observed, as found by Matos, Freitas, and Lo Monaco (2010) in cultivated HSSF-CWs for treating swine wastewater, and an increase was observed by Prata, Matos, Cecon, Monaco, and Pimenta (2013) in HSSF-CWs used to treat SW.
The correlation between d in HSSF-CWs and the n parameter in MFOKM, obtained for the COD and BOD variables, is presented in Table 4.The Pearson coefficient presents values ranging from -1 (negative correlation) to 1 (positive correlation), while values close to 0 indicate low correlation between the variables.Thus, a higher correlation between d and n was verified for COD values than for BOD values, where a smaller correlation was obtained for values of n from the soluble fraction of COD.Based on these results, it may be considered that the d variable and the n parameter are inversely correlated, and that the dispersion conditions of the environment influence the values of n, confirming the suspicions of Brasil et al. (2007).
Suspecting that the n value might also be related to the physical and biological processes, which occur in great intensity at the beginning of the HSSF-CWs, it can be interesting to observe the correlation between n and the decay intensity of the C/C 0 relationship, which happens in the first fifth of the HSSF-CWs length, a variable denominated Δ(C/C 0 ) 1/5 .The n values were then plotted in function of the Δ(C/C 0 ) 1/5 values, obtained in different HSSF-CWs for total and soluble COD and BOD, with the obtained curve presented in Figure 3.
As can be observed, n is strongly correlated to the decay in organic matter which occurs at the beginning of the HSSF-CWs, independent of the nature of the organic matter (BOD or COD, total or soluble).This correlation is shown both graphically as well as by the Pearson coefficient, which resulted in values of -0.985, indicating the existence of a strong inversely proportional correlation between the studied variables.In other words, the greater the variation magnitude in the concentrations of organic matter, at the beginning of HSSF-CWs, lower is the n value.
In Figure 4 the n x Δ(C/C 0 ) 1/5 graph is presented, plotted using data obtained in this work, along with data from Brasil et al. (2007) and Chagas et al. (2011) obtained in HSSF-CWs with L W -1 = 24.Even without obtaining an elevated R² value, maintaining of the inverse behaviour of n in relation to Δ(C/C 0 ) 1/5 is observed despite the diverse set of data, obtained in distinctive operational conditions and also in different HSSF-CWs.Based in the results obtained, an inverse correlation was verified both in relation to the efficiency of organic matter removal at the beginning of the HSSF-CWs as well as the dispersion in this medium with the HRT power parameter of the MFOKM adjusted to the data obtained in the system.

Conclusion
The different geometric configurations did not generate great differences with respect to overall removal of organic matter and dispersion among the HSSF-CWs evaluated.All HSSF-CWs demonstrated disperse behaviour according to the hydraulic performance.For this reason the equations which consider that the piston flow regime occurs during drainage in HSSF-CWs, as in the conventional first order kinetic, might be considered less adequate to describe the behaviour of these reactors.The power parameters of HRT in the MFOKM adjusted to the data obtained in the system (n), and the relative decay in concentration of organic matter in the first fifth of the HSSF-CWs lengths (Δ(C/C 0 ) 1/5 ), retained an inverse correlation with the dispersion coefficient (d).

Figure 1 .
Figure 1.HSSF-CWs sketch.a) dimensions and positioning details of the sample collection positions, with one of detailed individually.b) representation in a longitudinal section AA' with details of the hydraulic systems: influent inlet and effluent drain.c) sketch of the SW distribution and sample collection tube.

Figure 3 .
Figure 3. Graph and linear equation adjusted for n obtained after the adjustment of MFOKM of the data in function of C/C0 variation values in the first fifth of the HSSF-CWs length (Δ(C/C0)1/5).

Table 1 .
Mean characteristics and standard deviation of the raw sewage used in the experiment.

Table 2 .
Constructive and operational characteristics of the HSSF-CWs.

Table 3 .
Best-fitting model parameter for the Modified First Order Kinetic Model (MFOKM) obtained with the concentrations corrected for water losses by evaporation, together with the associated Coefficients of Determination (R 2 ).

Table 4 .
Pearson's linear correlation between the values of dispersion number (d) and the n parameter of the MFOKM adjusted to the data obtained in HSSF-CWs with different geometric configurations.