Sensitivity Analyses for Modeling Evolving Reactivity of Granular Iron for the Treatment of Trichloroethylene

To better predict long-term performance of a remediation system, parameters of a numerical model should be constrained with care by calibrating with reliable experimental data. This study conducted sensitivity analyses for model parameters, which were shown to represent reasonably well the observed geochemical behaviors for the column experiments that evaluated evolving reactivity of granular iron for the treatment of trichloroethylene (TCE) resulting from precipitation of secondary minerals. The particular model parameters tested include iron corrosion rate, aragonite and Fe2(OH)2CO3 precipitation rates, and proportionality constants for each mineral. For sensitivity analyses, a specific parameter was systematically changed, while other parameters were fixed at the values for the base case. The ranges of parameters tested were determined based on the previous modeling study. The results showed that the most important and sensitive model parameters were secondary mineral precipitation rates. Also, not only absolute precipitation rate for each mineral but also relative precipitation rates among different minerals were important for system performance. With help of sensitivity analysis, the numerical model can be used as a predictive tool for designing an iron permeable reactive barrier (PRB) and can provide implications for the long-term changes in reactivity and permeability of the system.


Introduction
Permeable reactive barriers (PRBs) using granular iron as a reactive material have been successfully applied for the removal of various contaminants in groundwater [1][2][3], and they are now accepted as one of the standard remediation options for contaminated groundwater [4,5].
While field demonstrations have shown that PRBs can be effective in the long term [1,6], concerns remain how long reactivity of iron and permeability of the systems can last under various geochemical conditions.Studies have shown that major geochemical components, such as calcium and carbonate, and oxidants, such as nitrate and chromate, can negatively affect the reactivity and permeability of iron PRBs [7][8][9].
Laboratory [10,11], field [1,6,12], and reactive transport modeling [13][14][15] experiments have been conducted to evaluate the reactivity and permeability of iron PRBs.In particular, reactive transport modeling can be used to assess various geochemical conditions without having extensive labor and cost for changing various conditions for laboratory and field tests.Given reasonable model parameter sets, the established model can be used to predict future performance of the systems.
Jeen et al. [16], in particular, showed that reactive transport modeling can be effective to demonstrate the geochemical patterns observed over time in a column experiment.A multi-component reactive transport model was modified to reflect the evolving reactivity of iron caused by precipitation of secondary carbonate minerals on the iron surfaces.The key consideration of the model was to update the reactive surface area of the iron, as geochemical conditions changed temporally and spatially over time.The modified model has since been successfully used to demonstrate the observed behaviors of geochemistry over time for laboratory experiments [8,[17][18][19] and field tests [20].
One of the remaining concerns about simulating long-term performance of iron PRBs using a model, as an example of Jeen et al. [16], is the predictive capability of the model with regard to the sensitivity of model parameters.By evaluating more sensitive model parameters for predictive simulation, one can spend more effort to get reliable model parameter values to improve predictability.Thus, the purposes of this study were to evaluate the sensitivity of model parameters for modeling evolving reactivity of iron, to identify more important parameters to predict the long-term performance of iron PRBs, and thus to provide insights for evolution of the systems with regard to changes in reactivity and permeability.This study conducted sensitivity analysis for model parameters, which were shown to be important to represent the geochemical evolution of a laboratory column experiment.

Model
The numerical model MIN3P [21] was developed for multi-component reactive transport simulations in variably saturated porous media.The governing equations were discretized using a locally mass conservative finite volume method and were linearized by a global implicit solution method [22,23].The formulation of the model was based on a partial equilibrium approach [22,24,25], allowing inclusion of kinetic limitations, which are essential for adequate description of groundwater treatment by PRB systems.
Previous studies showed that the reaction of granular iron and groundwater causes precipitation of secondary minerals, and the precipitation front migrates as the reactivity of iron is decreasing because of the deposition of secondary minerals on iron surfaces [26,27].In order to adequately reproduce the changes in iron reactivity, it is key to incorporate the evolving reactivity of the iron into the kinetic formulations in MIN3P.Inclusion of the evolving reactivity of iron was realized by updating the reactive surface area of the iron based on the laboratory-derived relationship between secondary mineral volume fraction and reactive surface area of iron [16].
In MIN3P, the rate expression for dissolution-precipitation of minerals is described as follows: where k eff,i is the effective rate constant for the dissolution of mineral phase i (mol L −1 H 2 O s −1 ), S i is the reactive surface area of mineral phase i (m 2 L −1 bulk), and IAP m i and K m i are the ion activity product and the corresponding equilibrium constant for mineral dissolution-precipitation reactions, respectively.The two important secondary minerals identified by X-ray diffraction (XRD) [28] and tested for the sensitivity analyses in this study are aragonite and Fe 2 (OH) 2 CO 3 .
The extent to which the accumulation of an individual mineral phase affects the reactivity of iron is represented by: where S(x, t) is the reactive surface area of the iron at a specific location along the flow path and time (m 2 iron L −1 bulk), S 0 is the initial reactive surface area of iron (m 2 iron L −1 bulk), α i is the proportionality constant for mineral phase i, and ϕ i (x, t) is the volume fraction of the mineral phase i at a specific location and time (-).The proportionality constant for each mineral phase represents the extent to which a particular mineral phase contributes to the reactivity loss of iron.
The reaction rate for trichloroethylene (TCE) degradation is represented by: where k SA-TCE-Fe 0 is the rate constant of TCE normalized to iron surface area (mol m -2 iron s −1 ), [TCE] is the concentration of TCE (mol L −1 H 2 O), and K 1/2 is the half-saturation constant (mol L −1 H 2 O), which corresponds to the TCE concentration at half-maximum transformation rate.It was considered that about 5% of TCE is degraded by hydrogenolysis to cis-1,2-dichloroethylene (cis-1,2-DCE) and the remaining 95% is degraded to ethene via β-elimination [16].The rate expression for iron corrosion by water is represented by: where k SA-H 2 O-Fe 0 is the rate constant of iron corrosion normalized to iron surface area (mol m −2 iron s −1 ), IAP H 2 O-Fe 0 is the ion activity product, and K H 2 O-Fe 0 is the equilibrium constant.In MIN3P, porosity and hydraulic conductivity are updated to reflect the changes during the simulation as secondary minerals precipitate and iron dissolves.The hydraulic conductivity was updated based on a normalized version of the Kozeny-Carmen relationship.The details of the model description and model validation processes can be found in Jeen et al. [16], and description of the column test, to which the model was applied, can be found in Jeen et al. [28].The sensitivity analyses in this study were performed for column D of Jeen et al. [16,28].

Model Parameters
Based on the observed geochemical changes in the column, the following 11 aqueous components were included for the expression of all chemical reactions: Ca 2+ , Cl − , CO 3 2− , Fe 2+ , Fe 3+ , H + , H 2 (aq), TCE, cis-1,2-DCE, vinyl chloride (VC), ethene, and H 2 O.A total of 12 aqueous complexes were also included for appropriate determination of mineral solubilities.The equilibrium constants were taken from the databases of MINTEQA2 [29] and WATEQ4F [30].
To represent the column system, one-dimensional discretization was used, with a spatial discretization interval of 0.01 m for the 0.5 m long column, giving a total of 51 grid points.The flow system was modeled as a fully saturated system with a second type (specified flux) boundary at the influent end and first type (specified head) boundary at the effluent end.The flow rate (1.38 × 10 −5 m s −1 ) at the influent of the column was used for the specified flux at the influent end, and zero hydraulic head was specified at the effluent end.Zero hydraulic head for the entire domain was specified as the initial condition for the flow simulation.
A diffusion coefficient of 1.5 × 10 −9 m 2 s −1 [31] and longitudinal dispersivity of 9.9 × 10 −4 m were used as the transport parameters.Boundary conditions for reactive transport were the third type (specified mass flux) at the influent end and the second type (free exit) at the effluent end.The measured aqueous concentrations and pH of the source water for each column were used as the influent chemical compositions.The initial condition for reactive transport, several orders of magnitude lower concentrations for each constituent, relative to the expected concentrations for the simulation, were specified in the entire domain.The detailed input parameters for the column test are summarized in Table 1.

Sensitivity Analysis
For constraining the model parameters in Jeen et al. [16], adjustment of parameters was conducted to fit overall profiles over time, not to fit a single profile at a specific time.As a result, each of the model's parameters had some ranges between different columns.The iron corrosion rate, aragonite and Fe 2 (OH) 2 CO 3 precipitation rates, and proportionality constants for each mineral were those with some ranges, and thus were parameters tested with sensitivity analyses in this study.
One set of model parameters, which provided the best fit to the data of column D, was selected as the base case for the sensitivity analyses.In Jeen et al. [31], column D received the highest concentration of aqueous calcium carbonate, showed the highest measured iron corrosion rate, and thus showed the fastest passivation of the iron over time [31].Therefore, any changes in the model parameters could be reflected most sensitively to the modeling results.
For sensitivity analyses, a specific parameter was systematically changed, while other parameters were fixed at the values for the base case.The parameter values for the base case and the range tested for each parameter are shown in Table 2.
Table 2.The parameter values for the base case and range tested for each parameter.

Parameter Base Case Sensitivity Analysis
Initial reactive surface area (m 2 iron L −1 bulk) 8.37 × 10 3 Fixed

Iron Corrosion Rate
The range selected for the initial iron corrosion rate was +/−20%, +/−50% from the base case.The fitted iron corrosion rates were generally very close to the measured maximum iron corrosion rates [16].Therefore, it was considered that the selected range for sensitivity analyses covers the possible changes induced by differences in the initial iron corrosion rate.Figure 1 shows the effect of the iron corrosion rate on the TCE profiles at different times.The deviations from the base case increased with time, approaching a maximum on Day 164.However, the effect leveled off at later times (Figure 1f) because of the exponential nature of the reactivity loss (Equation ( 2)).
Water 2018, 10, x FOR PEER REVIEW 6 of 17 With regard to mineral precipitation, when some amount of aragonite had accumulated, further increase in precipitation did not affect the reactivity to the same extent as the previous increase.The shape of accumulating aragonite did not change over time for different values of iron corrosion rate, but the amount of aragonite was proportional to the values of the iron corrosion rate (Figure 2).

Aragonite Precipitation Rate
The range of the fitted CaCO3(s)(aragonite) precipitation rate between different columns was within 1 order of magnitude [16].However, the range for sensitivity analyses was selected to be over 2 orders of magnitude because the rate constant for CaCO3(s) reported in the literature [32] was about 12 orders of magnitude lower than the range of the fitted aragonite precipitation rate.
Figures 3 and 4 show the effect of aragonite precipitation rate on the TCE profiles and aragonite volume fraction at different times, respectively.The effect of the aragonite precipitation rate was prominent over the range at which that the sensitivity analyses were conducted, as seen in the differences in the curvature of the TCE profiles that evolved over time (Figure 3).TCE was degraded much faster at a slower aragonite precipitation rate because of a slower increase in aragonite accumulation and thus a slower decrease in reactivity of the iron.
Aragonite precipitated further into the column at a slower precipitation rate, while aragonite precipitated closer to the influent end at a faster precipitation rate (Figure 4).However, eventually more aragonite was accumulated throughout the column at a faster aragonite precipitation rate

Aragonite Precipitation Rate
The range of the fitted CaCO 3 (s)(aragonite) precipitation rate between different columns was within 1 order of magnitude [16].However, the range for sensitivity analyses was selected to be over 2 orders of magnitude because the rate constant for CaCO 3 (s) reported in the literature [32] was about 1−2 orders of magnitude lower than the range of the fitted aragonite precipitation rate.
Figures 3 and 4 show the effect of aragonite precipitation rate on the TCE profiles and aragonite volume fraction at different times, respectively.The effect of the aragonite precipitation rate was prominent over the range at which that the sensitivity analyses were conducted, as seen in the differences in the curvature of the TCE profiles that evolved over time (Figure 3).TCE was degraded much faster at a slower aragonite precipitation rate because of a slower increase in aragonite accumulation and thus a slower decrease in reactivity of the iron.
Aragonite precipitated further into the column at a slower precipitation rate, while aragonite precipitated closer to the influent end at a faster precipitation rate (Figure 4).However, eventually more aragonite was accumulated throughout the column at a faster aragonite precipitation rate (Figure 4e,f).Accumulation of mineral was dependent not only on the absolute value but also on the relative precipitation rate among other mineral phases.Therefore, at a slower aragonite precipitation rate, there was greater accumulation of Fe 2 (OH) 2 CO 3 near the influent end (not shown).Because the molar volume of Fe 2 (OH) 2 CO 3 is greater than that of aragonite, more porosity loss occurred in that region.
Water 2018, 10, x FOR PEER REVIEW 8 of 17 (Figure 4e,f).Accumulation of mineral was dependent not only on the absolute value but also on the relative precipitation rate among other mineral phases.Therefore, at a slower aragonite precipitation rate, there was greater accumulation of Fe2(OH)2CO3 near the influent end (not shown).Because the molar volume of Fe2(OH)2CO3 is greater than that of aragonite, more porosity loss occurred in that region.

Fe2(OH)2CO3 Precipitation Rate
The range for sensitivity analyses of Fe2(OH)2CO3 precipitation rate was also selected to be over 2 orders of magnitude, based on the range of the fitted Fe2(OH)2CO3 precipitation rate between different columns [16].Generally, the effect of the Fe2(OH)2CO3 precipitation rate was similar to that of aragonite, but in the opposite direction and to a lesser extent (Figures 5 and 6).That is, at a slower Fe2(OH)2CO3 precipitation rate, greater reactivity loss occurred because of more precipitation of aragonite, resulting in further migration of the TCE profiles (Figure 5).Because the Fe2(OH)2CO3 precipitation rate indirectly affects the aragonite precipitation rate through a difference in relative precipitation rate, the extent to which Fe2(OH)2CO3 precipitation rate affects the TCE profiles was much less than that of aragonite (Figures 3 and 5).
For the accumulation of aragonite, at a slower Fe2(OH)2CO3 precipitation rate, more aragonite was accumulated near the influent end of the column (Figure 6).Similar to the case of aragonite

Fe 2 (OH) 2 CO 3 Precipitation Rate
The range for sensitivity analyses of Fe 2 (OH) 2 CO 3 precipitation rate was also selected to be over 2 orders of magnitude, based on the range of the fitted Fe 2 (OH) 2 CO 3 precipitation rate between different columns [16].Generally, the effect of the Fe 2 (OH) 2 CO 3 precipitation rate was similar to that of aragonite, but in the opposite direction and to a lesser extent (Figures 5 and 6).That is, at a slower Fe 2 (OH) 2 CO 3 precipitation rate, greater reactivity loss occurred because of more precipitation of aragonite, resulting in further migration of the TCE profiles (Figure 5).Because the Fe 2 (OH) 2 CO 3 precipitation rate indirectly affects the aragonite precipitation rate through a difference in relative precipitation rate, the extent to which Fe 2 (OH) 2 CO 3 precipitation rate affects the TCE profiles was much less than that of aragonite (Figures 3 and 5).
For the accumulation of aragonite, at a slower Fe 2 (OH) 2 CO 3 precipitation rate, more aragonite was accumulated near the influent end of the column (Figure 6).Similar to the case of aragonite precipitation, eventually more aragonite was accumulated throughout the column at a slower Fe 2 (OH) 2 CO 3 precipitation rate (Figure 6e,f).

Proportionality Constant for Aragonite
The selected range for sensitivity analysis of a1 (proportionality constant for aragonite) was +/−20%, +/−50%, which covers the maximum range of the fitted constant between different columns [16].The proportionality constant for each mineral phase affects the extent to which the reactivity is decreased due to the accumulation of that mineral.Therefore, at a higher value of 1  , more migration of the TCE profiles occurred, as seen in Figure 7.In contrast, aragonite was accumulated more at a lower value of 1  , because a slower decrease in reactivity provided more potential for aragonite to accumulate in the same period of time (Figure 8).

Proportionality Constant for Aragonite
The selected range for sensitivity analysis of a 1 (proportionality constant for aragonite) was +/−20%, +/−50%, which covers the maximum range of the fitted constant between different columns [16].The proportionality constant for each mineral phase affects the extent to which the reactivity is decreased due to the accumulation of that mineral.Therefore, at a higher value of α 1 , more migration of the TCE profiles occurred, as seen in Figure 7.In contrast, aragonite was accumulated more at a lower value of α 1 , because a slower decrease in reactivity provided more potential for aragonite to accumulate in the same period of time (Figure 8).

Proportionality Constant for Fe2(OH)2CO3
Because the absolute value of the fitted 2  (proportionality constant for Fe2(OH)2CO3) was low (0.0-5.0 for all columns; [16]), the range of +/−100% for sensitivity analysis of 2  does not necessarily mean a wide range of the parameter.Figures 9 and 10 show the effect of 2  on the TCE profiles and aragonite volume fraction at different times, respectively.As expected with the low value of 2  , the effect of this parameter was not significant.At a higher value of 2  , slightly more migration of the TCE profiles occurred over time (Figure 9).In contrast to the case of 1  , this parameter is applied to the accumulated volume fraction of Fe2(OH)2CO3 instead of that of aragonite.Because the relative precipitation rates for each mineral do not change, the shape and pattern of aragonite accumulation was similar to the case of 1  .There was greater aragonite accumulation over time at a lower value of 2  , but the difference was not significant, as expected by the low value of 2  (Figure 10).Because the absolute value of the fitted α 2 (proportionality constant for Fe 2 (OH) 2 CO 3 ) was low (0.0-5.0 for all columns; [16]), the range of +/−100% for sensitivity analysis of α 2 does not necessarily mean a wide range of the parameter.Figures 9 and 10 show the effect of α 2 on the TCE profiles and aragonite volume fraction at different times, respectively.As expected with the low value of α 2 , the effect of this parameter was not significant.At a higher value of α 2 , slightly more migration of the TCE profiles occurred over time (Figure 9).In contrast to the case of α 1 , this parameter is applied to the accumulated volume fraction of Fe 2 (OH) 2 CO 3 instead of that of aragonite.Because the relative precipitation rates for each mineral do not change, the shape and pattern of aragonite accumulation was similar to the case of α 1 .There was greater aragonite accumulation over time at a lower value of α 2 , but the difference was not significant, as expected by the low value of α 2 (Figure 10).

Conclusions
To reliably predict long-term performance of a remediation system, it is necessary to constrain the parameters of a numerical model with validation of real data.Sensitivity analysis demonstrates how well a particular parameter is constrained and thus provides a basis for evaluation of geochemical profile changes over time.This study conducted sensitivity analyses for model parameters, which were shown to be important and also have some ranges to fit laboratory column data.Based on the previous modeling study, iron corrosion rate, aragonite and Fe2(OH)2CO3 precipitation rates, and proportionality constants for each mineral were evaluated.
The iron corrosion rate affected TCE degradation and aragonite accumulation in proportion to the values of the rate.Secondary mineral precipitation rates were considered to be the most important parameters for predictive purposes.Interestingly, aragonite and Fe2(OH)2CO3 precipitation rates affected each other in opposite ways, indicating that not only absolute precipitation rate for each mineral but also relative precipitation rates among different minerals are important for system

Conclusions
To reliably predict long-term performance of a remediation system, it is necessary to constrain the parameters of a numerical model with validation of real data.Sensitivity analysis demonstrates how well a particular parameter is constrained and thus provides a basis for evaluation of geochemical profile changes over time.This study conducted sensitivity analyses for model parameters, which were shown to be important and also have some ranges to fit laboratory column data.Based on the previous modeling study, iron corrosion rate, aragonite and Fe 2 (OH) 2 CO 3 precipitation rates, and proportionality constants for each mineral were evaluated.
The iron corrosion rate affected TCE degradation and aragonite accumulation in proportion to the values of the rate.Secondary mineral precipitation rates were considered to be the most important parameters for predictive purposes.Interestingly, aragonite and Fe 2 (OH) 2 CO 3 precipitation rates affected each other in opposite ways, indicating that not only absolute precipitation rate for each mineral but also relative precipitation rates among different minerals are important for system performance.With higher values for the proportionality constants, the reactivity of iron for the treatment of TCE declined faster over time; however, it ended up with less accumulation of secondary minerals.This study implies that the reactivity and permeability of iron PRBs are not affected in the same direction.
Model parameters are often calibrated rather than directly measured because it is not possible to collect all necessary information with existing data.In this case, the parameters have some uncertainties, and thus model prediction should be interpreted with great care.The merit of model prediction is that it gives us implications with regard to future performance rather than precise particular numbers for certain scenarios.This study provides implications for long-term performance of granular iron PRBs which are affected by secondary mineral precipitation.It shows that the model can be a useful design tool if it is constrained with real data that are carefully generated by a laboratory test.

Water 2018 ,
10, x FOR PEER REVIEW 10 of 17precipitation, eventually more aragonite was accumulated throughout the column at a slower Fe2(OH)2CO3 precipitation rate (Figure6e,f).

Table 1 .
Input parameters used in the simulation.
The iron corrosion rate constant normalized to iron surface area; d The effective rate constant for the dissolution of mineral phase i; e The proportionality constant for aragonite; f The proportionality constant for Fe 2 (OH) 2 CO 3 (s).
Notes: a The trichloroethylene (TCE) rate constant normalized to iron surface area; b Half-saturation constant; c