Mechanistic modelling of Fe3+-EDDS mediated photo-fenton revisited: lumped radicals and sensitivity analysis

Chemical Engineering Journal 464 (2023) 142559 Available online 22 March 2023 1385-8947/© 2023 Elsevier B.V. All rights reserved. Mechanistic modelling of Fe 3+ -EDDS mediated photo-Fenton revisited: Lumped radicals and sensitivity analysis Kourosh Nasr Esfahani a , Montserrat P´ erez-Moya a , Mois´ es Graells a , Sara Miralles-Cuevas b , Alejandro Cabrera-Reina b , * a Universitat Polit´ ecnica de Catalunya, Av. Eduard Maristany, 16, 08019 Barcelona, Spain b Programa Institucional de Fomento a la Investigaci´ on, Desarrollo e Innovaci´ on (PIDi), Universidad Tecnol´ ogica Metropolitana, Av. Ignacio Valdivieso 2409, Santiago, Chile ARTICLE INFO Keywords: Photo-Fenton Photo-catalytic reaction Kinetic model Microcontaminants Fe 3+ -EDDS Lumped radicals ABSTRACT This study proposes a new kinetic model of the Fe 3+ -EDDS mediated photo-Fenton process at circumneutral pH when applied to microcontaminants (MCs) removal. The model is also able to work in the absence of H 2 O 2 , when EDDS •3− radicals generated from the lysis of the Fe 3+ -EDDS complex are responsible for MCs degradation, representing a new advance in this research field. A novel semi-empirical approach based on lumping radical species is adopted. This is done by including a factor in the model to describe the lower oxidation capacity of EDDS •3− radicals in front of hydroxyl radicals. Model calibration demonstrated the effectiveness of the semiempirical strategy to successfully predict the system behaviour in terms of sulfamethoxazole (SMX, model MC), Fe 3+ -EDDS and H 2 O 2 (when needed) concentration evolution during the process. Then, a global sensitivity analysis (GSA) was carried out to reduce the computation cost of the model, indicating that once the initial Fe 3+ - EDDS complex is oxidized to Fe 3+ -EDDS ox , its photo-activation to (Fe 3+ -EDDS ox )* is not actually reversible despite possible. Finally, the model was validated, showing that experimental data could be predicted properly, with NRMSE (Normalized Root Mean Square Errors) <0.08 for SMX and <0.15 for Fe 3+ -EDDS (normalized data). The oxidation capacity of the EDDS •3− radicals was estimated to be approximately 11% that of the hydroxyl radicals. 1. Introduction Emerging microcontaminants (MCs) are synthetic or natural substances, which are biorecalcitrants to conventional biological processes applied in wastewater treatment plants (WWTPs). MC term includes a large range of chemical compounds that present very different chemical natures hindering their removal [1]. Usually, they are emitted from point and nonpoint sources and end up in low concentrations from ng/L to µg/L in the water environment [2]. Pharmaceuticals, as an example of MCs, have obtained growing attention due to the increasing hazard posed by the presence of these chemicals and/or their transformation products in the natural environment that, even in very low concentrations, threatens human health [3]. Among them, Sulfamethoxazole (SMX) is a common antibiotic prescribed to treat infectious and respiration diseases, which is frequently detected in surface water ecosystems [4]. Since SMX contains a sulfonamide group that is present in many drugs, it is a representative drug for a broad range of antibiotics [5]. This antibiotic has been classified as an emerging contaminant, and has been found in water resources at concentration levels from ng/L to µg/L [6,7]. Thus, considering that SMX is an archetypal and extended pharmaceutical that presents frequent occurrence in different water bodies, it has been selected as a model pollutant to study the removal of MCs in WWTP secondary effluents by different processes such as membrane technologies [8], heterogeneous photocatalysis [9] or even the combination of different treatments [10], among many others. One of the most promising technologies for the degradation of MCs is advanced oxidation processes (AOPs), and among all AOPs, the photoFenton process has shown a high oxidation efficiency [11,12]. Several studies have demonstrated that the photo-Fenton process efficiently degrades a wide range of MCs in different water matrices at lab-scale Abbreviations: MCs, microcontaminants; SMX, sulfamethoxazole; WWTP, wastewater treatment plants; EDDS, ethylenediamine disuccinic acid; OM, organic matter; IC, inorganic carbon; DOC, dissolved organic carbon; TIC, total inorganic carbon; SA, sensitivity analysis; RMSE, root mean square error. * Corresponding author. E-mail address: [email protected] (A. Cabrera-Reina). Contents lists available at ScienceDirect Chemical Engineering Journal journal homepage: www.elsevier.com/locate/cej https://doi.org/10.1016/j.cej.2023.142559 Received 18 January 2023; Received in revised form 16 March 2023; Accepted 20 March 2023 Chemical Engineering Journal 464 (2023) 142559 2 [13], pilot-scale [14] and even large-scale [15]. However, the photo-Fenton process is pH-dependent due to the insolubility of the ferric aquo or hydroxy species at circumneutral pH, which limits its full-scale application due to the necessity of acidic pH. Consequently, the possibility of performing the photo-Fenton process at near-neutral pH has been studied using different strategies [16–19]. In this regard, polycarboxylic compounds that form a complex with iron, such as citrate, oxalate, maleate, and ethylenediamine disuccinic acid (EDDS) have been investigated to operate the photo-Fenton process at near-neutral pH [20–23]. Probably, the most studied chelating agent is EDDS, which is a biodegradable structural isomer of Ethylenediaminetetraacetic acid (EDTA) with efficient performance in the pH range 3–9 [24–26]. The reduction of Fe 3+ -EDDS to Fe 2+ -EDDS is a crucial step that governs the formation of hydroxyl radicals, mainly responsible for the degradation of the contaminants [27] but, additionally, the use of EDDS under UV radiation causes the production of organic radicals that directly or indirectly contribute to MCs removal. From an analytical and computational point of view, the development of robust mathematical models governing the dynamics and transport phenomena is a crucial step to improving the deep knowledge of the system and the process efficiency of almost any existing process. The availability of kinetic models allows for eliminating the dependency on the time and energy needed for intensive experimentations, which can be substituted by process simulation [28]. The kinetic modelling of the photo-Fenton process presents an important difficulty due to its complex reaction mechanism that has resulted in a lack of non-empirical models. Despite this, important efforts have been made in this field, especially by the research groups of Alfano from Consejo Nacional de Investigaciones Científicas y T´ ecnicas (Argentina) [29,30], G. Li Puma from Loughborough University (UK) [31,32] or J.A. S´ anchez-P´ erez from the University of Almería (Spain) [26,33]. In this context, while the kinetic modelling of the photo-Fenton process when operated at acidic pH has been typically investigated [34–37], the literature reveals low records on kinetic modelling at circumneutral pH. Within this frame of reference, only the mechanistic modelling of Fe 3+ -EDDS and Fe 3+ -NTA mediated photo-Fenton has been assessed. Soriano-Molina et al. (2018) [26] studied the modelling of the photo-Fenton process with Fe 3+ -EDDS, reporting higher photon absorption of iron complexed with EDDS compared to iron aquo complexes [26]. The model performance against different perturbations was also evaluated by analyzing its prediction capability when the treatment was applied to remove micropollutants from five different Spanish MWWTP secondary effluents [38]. Recently, Gualda-Alonso et al. (2022) [39] studied the mechanistic modelling of solar photo-Fenton mediated by Fe 3+ -NTA when applied to the removal of microcontaminants. Despite both models presenting excellent results, none of them implicitly included the effect of the organic radicals generated from the lysis of the iron chelate on process performance, limiting the model application only to the addition/presence of the oxidant agent (hydrogen peroxide). This work faces up the mechanistic modelling of the Fe 3+ -EDDS mediated photo-Fenton process, including the stand-alone option of oxidant agent absence (UVA/Fe 3+ -EDDS process). With this purpose in mind, SMX removal in synthetic MWWTP secondary effluent by UVA/ Fe 3+ -EDDS process under different UVA radiation levels was initially evaluated. Afterwards, the same procedure was done but applying the Fe 3+ -EDDS mediated photo-Fenton. A new semi-empirical mechanistic model was then proposed, and the experimental data was used for parameter estimation and model prediction capability evaluation. Finally, the impacts of fluctuations of the dynamic model parameters on the fitting performance were assessed through a sensitivity analysis reviewing the model structure. 2. Materials and methods 2.1. Reagents and chemicals SMX (>99%), selected as a model pollutant, was provided by Sigma Aldrich. Ferrous source (Fe 2 (SO 4 ) 3 ⋅H 2 O) and hydrogen peroxide (H 2 O 2 , 33% w/v) used for the photo-Fenton process were supplied by Panreac. EDDS was provided by Sigma-Aldrich. HPLC grade acetonitrile, methanol, glacial acetic acid (>99%, w/v), ortho-phenanthroline (99%, w/w) and sulfuric acid (96%, w/v) were also obtained from Panreac. All solutions were prepared in ultrapure water, produced with a Milli-Q® water purification system (PURELAB Option Q-7) with a specific resistance of 18.2 MΩcm and<5 µg/L of DOC. 2.2. Water matrix All the experiments were done in synthetic WWTP secondary effluent, which included: CaSO 4 ⋅2H 2 O (60 mg/L), MgSO 4 (60 mg/L), KCl (500 mg/L), (NH 4 ) 2 SO 4 (23.6 mg/L), K 2 HPO 4 (7.0 mg/L), NaHCO 3 (76.2 mg/L), beef extract (1.8 mg/L), peptone (2.7 mg/L), humic salts (4.2 mg/L), sodium lignin sulfonate (2.4 mg/L), sodium lauryl sulfate (0.9 mg/L), acacia gum powder (4.7 mg/L), and Arabic acid (5.0 mg/L). The total organic carbon (TOC) concentration of this recipe is ≈10 mg/L, and the total inorganic carbon (TIC) concentration is ≈15 mg/L, high enough to buffer the pH during the treatment so that it is possible to work at natural pH. 2.3. Experimental work plan The first experimental series was focused on the removal of SMX by the Fe 3+ -EDDS/UVA process, i.e., in the absence of H 2 O 2 . With this purpose in mind, several assays were done changing the power of the LED lights from 25% up to 100%. All the experiments, Table 1, were done at natural pH using 0.1 mM Fe 3+ -EDDS at 1:1 Fe 3+ :EDDS molar ratio. In this way, the contribution of the EDDS •3− radicals generated due to the breakage of Fe 3+ -EDDS could be determined. Then, the same experimental series was carried out but this time, adding the oxidant agent, i.e., using Fe 3+ -EDDS mediated photo-Fenton process. Once again, the experiments were done at natural pH using 0.1 mM Fe 3+ - EDDS at 1:1 Fe 3+ :EDDS molar ratio and 1.5 mM of H 2 O 2 . Thus, the contribution of the Fe 3+ -EDDS/UVA process on the Fe 3+ -EDDS mediated photo-Fenton process could be studied and modelled. 2.4. Experimental setup and experimental procedure The experimental series was done in a UVA-LED photoreactor (Fig. 1) consisting of two Pyrex glass tubes (5 cm internal diameter) connected to a recirculation tank. The system included 1200 LEDs (365 nm) equally distributed in 8 stripes and symmetrically placed around the tubes (4 stripes/tube). The stripes were connected to a computer through a datalogger so that the power of the LEDs in each stripe could be independently selected. The working volume was 8L, and water flow (20 L/min) was obtained by a centrifugal pump (Pan World NH-40PX). Initially, the pilot plant was loaded with 8L of the synthetic MWWTP secondary effluent, and SMX (1 mg/L) was added. After 10 min of recirculation, 0.1 mM Fe 3+ -EDDS was added, and after another 5 min of homogenization, the oxidant agent (when required) was added (1.5 Table 1 Experimental series. Assay ID SMX (mg/L) ±0.05 Fe 3+ -EDDS (mM) ±0.02 H 2 O 2 (mM) ±0.01 UV (W/m −2 ) A1 1.00 0.12 0 14 Fitting A2 1.00 0.12 0 42 Fitting A3 1.00 0.12 0 56 Fitting B1 1.00 0.12 1.50 14 Fitting B2 1.00 0.12 1.50 42 Fitting V1 1.00 0.12 0 28 Validation V2 1.00 0.12 1.50 28 Validation K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 3 mM). Immediately after, LED lights were switched on. Samples were taken every 2 min during the first 10 min, then every 5 min until 20 min and, from that moment, every 10 min until the end of the experiments (60 min). 2.5. Analytical determinations All samples related to the quantification of SMX and Fe 3+ -EDDS were quenched and filtered through 0.22 µm PTFE filters (Millipore) before injection in a JASCO Extrema HPLC equipped with a Supelco Fig. 1. Experimental setup. Fig. 2. Model scheme and semi-empirical approach. Table 2 Reaction rates of the proposed kinetic model. REACTION RATES r1=k1.[Fe3+−EDDS](I I+1)(1) r9=k9.[Fe3+−EDDSOX](I I+1)(9) r2=k2.[(Fe3+−EDDS)*](2) r10 =k10.[(Fe3+−EDDSOX)*](10) r3=k3.[(Fe3+−EDDS)*](3) r11 =k11.[(Fe3+−EDDSOX)*](11) r4=k4.[Fe3+−EDDS][H2O2](4) r12 =k12.[O.M][HO⋅](12) r5=k5.[Fe2+−EDDS](5) r13 =k13.[SMX][HO⋅](13) r6=k6.[Fe2+][H2O2](6) r14 =k14.[IC][HO⋅](14) r7=k7.[Fe3+−EDDS][H2O2](7) r15 =k15.[H2O2][HO⋅](15) r8=k8.[Fe3+−EDDS][HO⋅](8) Table 3 Mass balances of the proposed kinetic model. MASS BALANCES d[Fe3+−EDDS] dt = − r1+r2+r4−r7−r8 (16) d[(Fe3+−EDDS)*] dt =r1−r2−r3 (17) d[Fe2+] dt =r5−r6+r11 (18) d[Fe2+−EDDS] dt =r3−r4−r5+r7 (19) d[Fe3+−EDDSOX] dt =r8−r9+r10 (20) d[(Fe3+−EDDSOX)*] dt =r9−r10 −r11 (21) d[H2O2] dt = − r4−r6−r7−r15 (22) d[HO⋅] dt =r4+r6−r8−r12 −r13 −r14 −r15 +x⋅(r5+r11)(0<x<1)(23) d[OM] dt = − r12 (24) d[SMX] dt = − r13 (25) d[IC] dt = − r14 (26) K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 4 LICHROSPHER RP18-5 column and coupled to a JASCO DAD detector (model MD-4015). For SMX determination, an isocratic method was used in which the mobile phase A was 95% methanol and the mobile phase B was 5% formic acid (0.1%). SMX retention time was 4.7 min, and the resulting area was quantified at 267 nm. For Fe 3+ -EDDS determination, an isocratic method was used in which the mobile phase A was 5% methanol, and the mobile phase B was 95% buffer solution (15 mM sodium formate, 2 mM tetrabutylammonium hydrogen sulphate adjusted to pH 4 with acetic acid). The Fe 3+ -EDDS retention time was 6.8 min and the resulting area was quantified at 240 nm. H 2 O 2 was measured following DIN38402H15 based on titanium (IV) oxysulfate, while dissolved iron concentration was determined following ISO 6332 based on 1,10-phenantroline. Dissolved Organic Carbon (DOC) and total inorganic carbon (TIC) concentrations were measured in a Shimadzu TOC-VCN analyzer. Samples were filtered through 0.45 μ m Nylon filters (Milipore). 2.6. Kinetic model The kinetic model proposed in this work revisits the reaction scheme published by Soriano-Molina et al. (2018) [26] with the main objective of expanding its applicability by considering the contribution of the organic radicals yielded by the breakage of Fe 3+ -EDDS. However, the direct inclusion of the EDDS •3− radical into the mechanistic model inexorably involves the differentiation between the microcontaminant degradation reactions by hydroxyl radicals (HO • ) and by EDDS radicals (EDDS •3− ). This is almost unaffordable from a computational point of view for real applications. Actual MWWTP secondary effluents contain hundreds of microcontaminants; therefore, each new radical considered in the model would be translated into hundreds of new reactions and, consequently, into hundreds of additional model parameters. To overcome this critical problem, this work proposes a new semi-empirical approach consisting of lumping the concentrations of oxidizing radicals into the concentration of a single species, the hydroxyl radical. Hence, EDDS •3− radicals are assumed to be instantly converted into HO • radicals (Fig. 2). At the same time, however, its concentration is also assumed to be reduced by a factor (named “x” in the model) that represents the lesser oxidation capacity of the EDDS •3− radical in front of the HO • radical, a factor that will be later calibrated when fitting the model to the experimental data. In this way, despite the complexity of the real mechanism, the EDDS •3− radical can be easily introduced in the system using a single additional parameter (the reduction factor, 0<x <1). Regarding the model structure (Fig. 2), very briefly, it can be divided into the reaction mechanism related to radicals generation and the reactions related to radicals consumption. In the former (R1-R11), the initial Fe 3+ -EDDS evolves into Fe(OH) 3 and radicals (both HO • and EDDS •3− ) due to its interaction with H 2 O 2 , UVA radiation and the own radicals. In the latter, radicals can (efficiently) react with the model microcontaminant (R13) or can be (inefficiently) consumed by organic matter (OM) (R12), inorganic carbon (IC) (R14) or H 2 O 2 (R15), giving rise to oxidized organic matter (MX), oxidized inorganic carbon (ICX) and decomposed H 2 O 2 , respectively. OM includes the amount of total organic carbon provided by the chelating agent. Notice that R8, which is part of the radicals generation mechanism, also comprises the consumption of radicals. As explained above, the model considers the instant conversion of EDDS •3− into HO • by applying the reduction factor (semi-empirical approach). On the other hand, the dynamic model is based on reaction rates (Table 2) and mass balances to the most relevant components of the model (Table 3). Reaction rates are calculated following conventional equations in which a kinetic constant is multiplied by the concentration of the corresponding reagents participating in the reaction. In the reaction rates where UVA is involved, a saturation expression for Fig. 3. Experimental data (-●- SMX, -■- Fe 3+ -EDDS, -▴- H 2 O 2 ) for the I =14 W/m 2 experiment (empty symbols) and average experimental data for the I =28 W/m 2 , I =42 W/m 2 and I =56 W/m 2 experiments (closed symbols). Experimental conditions: [SMX] =1.00 ±0.05 mg/L; [Fe 3+ -EDDS] =0.12 ±0.02; (A) Fe 3+ -EDDS/ UVA assay [H 2 O 2 ] =0 mM and (B) photo-Fenton assay [H 2 O 2 ] =1.5 ±0.01 mM. Table 4 Estimation of the kinetic parameter values for the proposed model. Kinetic constants Initial value Fitted value (average) SD k1(min−1) 5.80⋅10 −2 5.30⋅10 −2 1.50⋅10 −2 k2(min−1) 2.50⋅10 −1 3.11 1.20⋅10 −2 k3(min−1) 1.71⋅10 +1 8.97⋅10 +2 2.44⋅10 +2 k4(mM−1.min−1) 1.90⋅10 +3 1.10⋅10 +4 2.70⋅10 +3 k5(min−1) 3.57⋅10 +1 3.45 2.06 k6(mM−1.min−1) 4.56 2.98⋅10 +2 3.37⋅10 +1 k7(mM−1.min−1) 4.00⋅10 −1 1.44⋅10 −1 4.10⋅10 −2 k8(mM−1.min−1) 5.60⋅10 +7 4.10⋅10 +9 8.21⋅10 +8 k9(min−1) 2.11⋅10 +1 1.17 4.36⋅10 −1 k10(min−1) 1.05⋅10 +1 5.61⋅10 −1 1.78⋅10 −3 k11(min−1) 3.31 3.40⋅10 +2 3.59 k12(mM−1.min−1) 1.76⋅10 +4 8.41⋅10 +4 2.30⋅10 +2 k13(mM−1.min−1) 4.80⋅10 +4 8.80⋅10 +9 1.99⋅10 +9 k14(mM−1.min−1) 5.10⋅10 +5 5.16⋅10 +5 1.42⋅10 +4 k15(mM−1.min−1) 1.60⋅10 +6 4.40⋅10 +8 5.38⋅10 +7 x 5.00⋅10 −1 1.15⋅10 −1 5.40⋅10 −2 K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 5 irradiance was included according to the experimental results (see Section 3.1). A similar strategy was proposed by Rivas et al. (2015) [40]. These equations should be ideally substituted in the future by more accurate expressions based on the volumetric rate of photon absorption; nevertheless, this deserves a whole scientific publication due to its high complexity. The concentration of OM and IC in the reaction rates is represented by DOC concentration and TIC concentration, respectively. The mass balances for each component can then be easily determined considering the generation and consumption rates of the reactions where they participate. Special attention requires the HO • mass balance. Since the conversion of EDDS •3− into HO • is instantaneous, the HO • mass balance also includes the EDDS •3− generation rates (R5 and R11) that are multiplied by the reduction factor (x). 2.7. Model parameter estimation The set of reaction rates given by Eq. 1–15 and the set of ordinary differential equations given by Eq.16–26 was implemented in MATLAB/ Simulink® (version R2022a) and solved numerically using ode15s solver with a time variable-step. The concentration profiles of SMX, Fe 3+ -EDDS and H 2 O 2 were used to estimate the model parameters (model calibration) that include 15 kinetic constants and the reduction factor. This was done by solving a nonlinear multivariate optimization problem to minimize the sum of the squared differences between the model predicted values and the corresponding experimental data available (Eq. (27). minZ =∑ m j=1 (∑ n i=1( [SMX](i,j)-[SMX](i,j)  [SMX](1,j))2 +( [H2O2](i,j)-[H2O2](i,j)  [H2O2](1,j))2 +( [Fe3+−EDDS](i,j)-[Fe3+−EDDS](i,j)  [Fe3+−EDDS](1,j))2 (27) Where the circumflex stands for the experimental data, the subscript i for each sample at specific reaction times and the subscript j represents the different experiments. This problem was solved using the nonlinear least-squares method with the Levenberg-Marquardt algorithm available in the Estimator Toolbox of MATLAB/Simulink®. The Root Mean Square Error (RMSE), as the standard deviation of the residuals (prediction errors), was calculated by Eq. (28) in order to test the model reliability and the goodness of fit. The coefficient of determination R 2 (Eq. (29) is also presented for illustrative purposes because it has been widely used, though R 2 in nonlinear models is not equal to the regression sum-ofsquares plus the residual sum-of-squares, as in the case of linear regression [41]. RMSEk= (∑ Nk i=1(yik-yik y1k)2/Nk) √ √ √ √(28) Rk2=1−(∑ N i=1 (yi-yi)2/∑ N i=1 (yi-¯y)2)(29) Fig. 4. Experimental data (symbols) and profiles predicted by the proposed kinetic model (lines). Experimental conditions: [SMX] =1.00 ±0.05 mg/L; [Fe 3+ - EDDS] =0.12 ±0.02 Mm; I=14 W/m 2 (25%); (A) Fe 3+ -EDDS/UVA assay [H 2 O 2 ] =0 mM and (B) photo-Fenton assay [H 2 O 2 ] =1.5 ±0.01 mM. Table 5 The goodness of fit (assays A1 and B1). Assay ID UV intensity NRMSE R 2 (W/m 2 ) SMX Fe 3+ -EDDS H 2 O 2 SMX Fe 3+ -EDDS H 2 O 2 A1 14 0.070 0.089 NA 0.91 0.97 NA B1 14 0.126 0.153 0.075 0.92 0.88 0.93 K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 6 where yi and yi correspond to the measured and simulated values at the given time, respectively, and ¯ y is the mean of the measured data. 2.8. Sensitivity analysis Parameter sensitivity analysis (SA) allows for examining a mathematical model and assessing the robustness of its output with respect to parameter uncertainty [42]. For a kinetic model, SA reveals the relative importance of the different steps included in the model. The selection of a SA approach among the various reported studies in the literature was discussed in a previous paper [34]. Hence, Global sensitivity analysis (GSA) is adopted since it is widely used for models with multiple correlated outputs. GSA assigns the outputs uncertainty to the uncertainty in each input factor over their entire range of interest. In many complex and nonlinear spatial phenomena, input factors usually show interactions and assessing the impact of a model input on a given output, being other factors constant, results inappropriate. Thus, SA is considered global when all the input factors are varied simultaneously, and the sensitivity is evaluated over the entire range of each input factor. Hence, SA quantifies the importance of model inputs and their interactions with respect to model output [43]. Among the assortment of GSA techniques, mostly statistical and analytical, the combined Latin Hypercube Sampling-Partial Rank Correlation Coefficient (LHS-PRCC) method [44] is used, following the analysis of previous work [34] to analyze the correlation between the model parameters and the output. A model with P parameters and a model output y can be written as a function of all the model parameter values so that y=f(x1,x2,⋯,xP). A distribution of N samples (i=1,2,⋯,N) for each model parameter xk (k=1,2,⋯,P), will produce a set of P sample vectors Xk (1 ×N) and an associated vector Y (1×N) of model output values. Hence, the partial rank correlation coefficient for the k-the parameter (RXkY) is given by [45]: RXkY=cov(Xk,Y)  Var(Xk)Var(Y) √=∑N i=1(xik −¯x)(yi−¯y)  ∑N i=1(xik −¯x)2(yi−¯y)2 √(30) Fig. 5. Profiles of the intermediates unobserved variables predicted by the proposed kinetic model. Experimental conditions: [SMX] =1.00 ±0.05 mg/L; [Fe 3+ - EDDS] =0.12 ±0.02 Mm; I =14 W/m 2 (25%): (A) Fe 3+ -EDDS/UVA assay [H 2 O 2 ] =0 mM and (B) photo-Fenton assay [H 2 O 2 ] =1.5 ±0.01 mM. Fig. 6. Sensitivity analysis and correlation coefficients R XY for all the proposed model parameters (absolute value): (A) SMX, (B) Fe 3+ -EDDS, (C) H 2 O 2 and (D) Overall. K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 7 The value of the partial rank correlation coefficient, ranging from −1 to 1, rates the extent to which an increase in the parameter value produces an increase (or decrease) in the model output [44]. Thus, it is commonly assumed that obtaining values in the range −0.05<RXY < 0.05 reveal a negligible significance of a parameter in the model output. A model producing multiple outputs may be managed by focusing on a reduced set of aggregated output values Y providing the most insight into the model analysis. In this case, the conformity of the model to the available experimental data is the criterion for selecting as the output values of interest the sum of the square errors produced for all the measured variables (SMX, Fe 3+ -EDDS and H 2 O 2 ). Thus, the output values for the sensitivity analysis and the calculation of the partial rank correlation coefficients are the following: Y1=∑ i( [SMX]i− [SMX]i)2(31) Y2=∑ j( [Fe3+−EDDS]j−[Fe3+−EDDS]j)2(32) Y3=∑ k( [H2O2]k− [H2O2]k)2(33) 3. Results and discussion 3.1. SMX removal by UVA/Fe 3+ -EDDS process and photo-Fenton mediated by Fe 3+ -EDDS SMX degradation by the UVA/Fe 3+ -EDDS process was initially studied. In this experimental series, the radiation level was varied from 14 W/m 2 (25%) to 56 W/m 2 (100%) to understand how this variable affects process performance (Fig. 3) and, afterwards, translate this behaviour into the model. Photo-saturation effect was observed. The SMX removal efficiency improved with irradiance level only up to a certain level, from which process performance remained almost constant. This is directly related to the Fe 3+ -EDDS evolution during the process. Since the oxidant species responsible for SMX degradation are generated from the Fe 3+ -EDDS breakage and equivalent Fe 3+ -EDDS concentration profiles were observed from 28 W/m 2 (50%), this also resulted in equivalent SMX removal curves. A similar trend was observed during the photo-Fenton treatment. Both Fe 3+ -EDDS concentration and H 2 O 2 consumption profiles showed once again almost equivalent behaviours from 28 W/m 2 , which resulted in extremely similar SMX removals from this irradiance level. 3.2. Model fitting and parameter estimation Model fitting is a non-convex optimization problem for which multiple optima may exist. Thus, different sets of parameter values (local optima) may be found depending on the starting point used by the search algorithm. Indeed, this may also lead to further discussion on global optimization, but this is beyond the scope of this work since the solution attained will finally be validated on experimental tests. In this work, the parameter values reported by Soriano-Molina et al. (2018) [26] were used as initial guesses (Table 4). Regarding parameter × (eq.23) since it is bounded between 0 and 1, the initial value was set to 0.5. The model was fit to each of the experimental data sets using the MATLAB/Simulink® R2022a Optimization Toolbox™ (lsqnonlin). The adjustment of the model to all experimental data sets (Table 1) always resulted in fitting with NRMSE>0.16 for all the cases. Model fitting to each experimental data set produced a set of parameter values that were averaged; these values and the corresponding standard deviation are given in Table 4. In the reactions related to radicals generation, the differences found between the initial guesses and the fitted values of the model parameters were lower than one order of magnitude. Probably, these differences are Fig. 7. Sensitivity analysis and correlation coefficients R XY for all parameters of the model without H 2 O 2 (Fe 3+ -EDDS/UVA assays): (A) SMX, (B) Fe 3+ -EDDS and (C) Overall. Table 6 Estimation of the kinetic parameters for the reduced model. Kinetic constants Initial value Fitted values (average) SD k1(min−1) 5.30⋅10 −2 8.49⋅10 −2 1.15⋅10 −1 k2(min−1) 3.11 3.12 1.58⋅10 −2 k3(min−1) 8.97⋅10 +2 7.61⋅10 +3 4.93⋅10 +3 k4(mM−1.min−1) 1.10⋅10 +4 1.23⋅10 +4 1.71⋅10 +3 k5(min−1) 3.45 2.32 2.38 k6(mM−1.min−1) 2.98⋅10 +2 2.53⋅10 +2 2.81⋅10 +1 k7(mM−1.min−1) 1.44⋅10 −1 1.41⋅10 −1 6.91⋅10 −2 k8(mM−1.min−1) 4.10⋅10 +9 5.11⋅10 +9 1.11⋅10 +9 k9(min−1) 1.17 1.25 3.09⋅10 −1 k11(min−1) 3.40⋅10 +2 3.54⋅10 +2 2.38⋅10 +1 k12(mM−1.min−1) 8.41⋅10 +4 8.43⋅10 +4 7.55⋅10 +1 k13(mM−1.min−1) 8.80⋅10 +9 1.04⋅10 +10 1.22⋅10 +9 k14(mM−1.min−1) 5.10⋅10 +5 4.20⋅10 +5 2.26⋅10 +5 k15(mM−1.min−1) 4.40⋅10 +8 3.95⋅10 +8 1.17⋅10 +8 x 1.15⋅10 −1 1.09⋅10 −1 2.02⋅10 −2 K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 8 related to the UV radiation source (solar light vs LED light) and how the catalyst-UV radiation interaction is mathematically considered in the model. On the other hand, the differences were higher in the reactions related to radicals consumption, which is indeed coherent because (i) the selected model MC was not the same in both works (acetamiprid vs SMX, being the former extremely recalcitrant) and (ii) the model proposed by Soriano-Molina does not consider the contribution of EDDS •3− radical to MCs degradation. In this sense, model fitting suggests that one EDDS •3− radical could be considered as ≈0.1 HO • radical. The SD values were found, in general, high but in the range of the ones obtained in similar modelling works [33]. Given the set of averaged model parameter values fit to the experimental data (Table 1), simulations were performed to check and illustrate the fitting. Fig. 4 shows the experimental data and predicted profiles for the same experimental conditions, with and without the addition of hydrogen peroxide (assays A1, B1): The goodness of fit (assays A1 and B1) is tabulated in Table 5. Despite their variability, the average values estimated for the model parameters (Table 4) allow producing practical model responses matching the experimental results. Fig. 4 compares the experimental data and the SMX, Fe 3+ -EDDS and H 2 O 2 concentration profiles obtained using such averaged parameter values. Fig. 4 compares two specific cases with (Assay B1) and without (Assay A1) H 2 O 2 addition and shows the capability of the model to reproduce the measured behavior and explain the SMX, Fe 3+ -EDDS and H 2 O 2 concentration profiles in these two extreme situations. Additional simulations are carried out to check the consistency of the model and the capacity to estimate the behavior of the intermediate unobserved variables, Fig. 5. Once the proposed model has revealed its practicality and capability to explain the process in a range of conditions beyond that addressed by the model by Molina et al. (2018) [26], it should be further examined in regard of parameter uncertainty and the need for the different mechanisms proposed in the model, which may reveal the chances for the simplification of the model. This is next performed following the scheme reported by Nasr Esfahani et al. (2022) [34]. 3.3. Global sensitivity analysis The complexity of the proposed model is next analyzed by means of Global Sensitivity Analysis (GSA). Model parameters are assumed independent, and a set of 500 samples for each model parameter is taken to apply the Latin Hypercube Sampling-Partial Ranking Correlation Coefficient (LHS-PRCC) analysis. The outputs considered for this analysis is the sum of the square error between measured and estimated values (SMX, Fe 3+ -EDDS and H 2 O 2 ), which focuses the analysis on the accuracy and fidelity of the model to the measured variables. Fig. 8. Experimental data (symbols) and profiles predicted by the reduced kinetic model (lines) after GSA. Experimental conditions: [SMX] =1.00 ±0.05 mg/L; [Fe 3+ -EDDS] =0.12 ±0.02 Mm; Fe 3+ -EDDS/UVA assays [H 2 O 2 ] =0 mM at (A) 14 W/m 2 ; (B) 42 W/m 2 and (C) 58 W/m 2 ; photo-Fenton assays [H 2 O 2 ] =1.5 ±0.01 mM at (D) 14 W/ m 2 and (E) 58 W/m 2 . K. Nasr Esfahani et al. Chemical Engineering Journal 464 (2023) 142559 9 The analysis is performed at two complementary levels. On the one hand, a partial GSA is performed using assays B1 and B2 with hydrogen peroxide, Fig. 6. On the other hand, is admitted that the presence of hydrogen peroxide may favor some reactions (R4, R6, R7 or R15), while masking some other (R5 and R11), which may result in mistakenly concluding the scarce relative importance of the related kinetic constants (model parameters). Hence, the complementary analysis, Fig. 7, is carried out using assays A1, A2 and A3, for which the absence of hydrogen peroxide will reveal the true significance of these mechanisms. It is worth noting that Fig. 7 shows non-zero values for the correlation coefficients of the parameters involved in the reactions that cannot occur in the absence of hydrogen peroxide (k4, k6, k7, k15). The values are very low, and they can be neglected, but the fact that the values are not strictly zero reveals the limitations of the LHS. While practical in terms of computational cost, LHS produces random samples for which the covariance between a set of random outputs and a set of noncorrelated inputs cannot be expected to be strictly zero. A particular one-at-a-time sampling for k15 confirmed the absolutely null influence of this parameter in the fitting of assays B1 and B2, for which there is no hydrogen peroxide. Fig. 6 and Fig. 7 show the absolute value of the averaged correlation coefficients found for all model parameters. The values are displayed in order of its importance for the fitting of the estimated values to the measured data available for SMX, Fe 3+ -EDDS and H 2 O 2 . Particularly, Fig. 6 reveals that correlation coefficients for k10,k5, and k6 are below the threshold |RXY| ≥0.05 for all three outputs SMX, Fe 3+ -EDDS and H 2 O 2 , which indicates their scarce relevance to the capability of the model to fit the experimental data. In addition, Fig. 7 reveals that correlation coefficients for k10,k12 , and k14 are below the threshold |RXY| ≥0.05 for all outputs SMX and Fe 3+ -EDDS, which indicates their scarce relevance to the capability of the model to fit the experimental data. The negligible effect of k10 in both scenarios suggests that the reverse reaction 10 is negligible for the global process. Consequently, reaction 10 is eliminated, assuming that once the initial Fe 3+ -EDDS complex is oxidized to Fe 3+ -EDDS ox , its photo-activation to (Fe 3+ -EDDS ox )* is not reversible. Certainly, this analysis and this assumption are acceptable only in the presence of light. The resulting simplified model was fit again to each of the experimental data sets as previously reported using as initial guesses the obtained values during the previous model calibration. The new average values obtained are reported in Table 6, along with the corresponding standard deviation. This interesting result shows that different model fittings could explain the same experimental data. This is clearly connected to the prior discussion on the multiple local optima that can be encountered by the optimization algorithm and reveals that while a model can be correctly adjusted for practical purposes, further mathematical research beyond the scope of this work is required to determine the best number of model parameters and their optimal value. Simulations were performed to check and illustrate the fitting after the GSA, predicted profiles and experimental data are overlapped for each experimental condition (Fig. 8). The simulations by the simplified model showed good fitting of the experimental data for both Fe 3+ - EDDS/UVA process and Fe 3+ -EDDS mediated photo-Fenton process. Fig. 9. Experimental data (symbols) and profiles predicted by the reduced kinetic model (lines) after GSA for the validation experiments. Experimental conditions: [SMX] =1.00 ±0.05 mg/L; [Fe 3+ -EDDS] =0.12 ±0.02 Mm; I =28 W/m 2 (50%); (A) Fe 3+ -EDDS/UVA assay [H 2 O 2 ] =0 mM with (B) the corresponding residuals and (C) photo-Fenton assay [H 2 O 2 ] =1.5 ±0.01 mM with (D) the corresponding residuals. K. Nasr Esfahani et al.