Transport and Retention of Sulfidized Silver Nanoparticles in Porous Media: The Role of Air‐Water Interfaces, Flow Velocity, and Natural Organic Matter

The sulfidation and aging of silver nanoparticles (Ag‐NPs) with natural organic matter (NOM) are major transformation processes along their pathway in wastewater treatment plants and surface waters. Although soils appear to be a sink for disposed Ag‐NPs, the impact of variable saturation on the transport and retention behavior in porous media is still not fully understood. We studied the behavior of sulfidized silver nanoparticles (S‐Ag‐NPs, 1 mg L−1) in saturated and unsaturated sand columns regarding the effects of (i) the presence of NOM (5 mg L−1) in the aquatic phase on retention, transport, and remobilization of S‐Ag‐NPs and (ii) the distribution and quantity of air‐water and solid‐water interfaces for different flow velocities determined via X‐ray microtomography (X‐ray μCT). Unsaturated transport experiments were conducted under controlled conditions with unit gradients in water potential and constant water content along the flow direction for each applied flux. It was shown that (i) NOM in S‐Ag‐NP dispersion highly increased the NP‐mobility; (ii) differences between saturated and unsaturated transport were increasing with decreasing flux and, consequently, decreasing water contents; (iii) both, solid‐water and air‐water interfaces were involved in retention of S‐Ag‐NPs aged by NOM. Using numerical model simulations and X‐ray μCT of flow experiments, the breakthrough of Ag‐NP could be explained by a disproportional increase in air‐water interfaces and an increasing attachment efficiency with decreasing water content and flow velocity.


Introduction
Silver nanoparticles (Ag-NPs) are among the most commonly used engineered NPs. Their wide application in detergents, clothing, cosmetics, and medical care products leads to a wide environmental spreading (Gottschalk et al., 2013). It is expected that a large part of engineered nanoparticles are released into aquatic and soil systems through wastewater treatment plants (Bundschuh et al., 2018). Concentrations of Ag in environmental compartments have been determined to be in the range of ng L −1 for water and mg kg −1 for soil and sediments (Blaser et al., 2008;Gottschalk et al., 2009;McGillicuddy et al., 2017). Along their pathway in sewer systems, wastewater treatment plants, and surface waters Ag-NP are sulfidized (S-Ag-NP) and coated by natural organic matter (NOM) (Brunetti et al., 2015;Dale et al., 2013;Kim et al., 2010;Lowry et al., 2012;Thalmann et al., 2016). NOM occurs in almost all aquatic systems in a wide concentration range of 0.1 to 60 mg L −1 (Frimmel & Abbt-Braun, 1999;Philippe & Schaumann, 2014). One of the most common classes of NOM is humic substances, which consist of organic macromolecules with various molecular weights that are negatively charged under natural conditions. For different NPs, it has been shown that aging in the presence of humic substances can modify their surface properties like wettability and surface charge, thereby in most cases reducing their tendency to aggregate (Akaighe et al., 2013;Baalousha et al., 2008;Christian et al., 2008). The transformation of Ag-NPs in the environment can significantly change their stability, mobility, and toxicity (Baalousha et al., 2016;Grillo et al., 2015).
The transport of engineered inorganic NPs through porous media is controlled by several key factors (El Badawy et al., 2013): the transport velocity and tortuosity, grain size and grain surface roughness Sagee et al., 2012), solution chemistry of the mobile phase (Makselon et al., 2017), geochemistry of the porous substrate, and the concentration, stability, surface charge, and aggregation behavior of NPs themselves. In quartz sand and natural soils, the presence of dissolved organic matter enhanced the mobility of stabilized Ag-NP (Sagee et al., 2012;Yecheskel et al., 2016) while low pH, higher ionic strength and divalent cations promoted NP aggregation and retention (Akaighe et al., 2013;Baalousha et al., 2008). In saturated, pure quartz sand with a negative surface charge Ag-NPs stabilized colloidally often show a high mobility whereas in natural soils, the decreasing grain size, the presence of clay minerals and positively charged iron and aluminum oxides causes enhanced NP retention (Cornelis et al., 2013;El Badawy et al., 2013;). An increasing flow rate and a decreasing surface area of potential collectors reduces the probability of NP attachment (McDowell-Boyer et al., 1986).
In the unsaturated zone, the presence of air as an additional phase enhances the complexity of the porous media influencing also particle transport (Flury & Aramrak, 2017). Depending on the water content, different pore sizes are involved in water flow, which control both, the transport regime (velocity and dispersion), and the available solid-water interface (SWI) as a potential collector for NPs (Knappenberger et al., 2014). Further, the tortuosity is increased by decreasing water contents, which increases the length of the flow paths and the probability for particle attachment. The air-water interfaces (AWI) provide an additional potential attachment surface, which is negatively charged (ζ-potential of −20 to −40 mV) at a moderate pH of 5.5 to 6.0 and ionic strength of 0 to 0.1 mmol L −1 NaCl (Leroy et al., 2012;Takahashi, 2005). Besides the weak electrostatic attraction for positively charged particles, the superhydrophobic AWI attracts hydrophobic surfaces of the particles (Ducker et al., 1994). For citrate stabilized Ag-NPs (Kumahor, Hron, et al., 2015) and soil-aged Ag-NPs (Kumahor et al., 2016) it was suggested that the AWI provides a reversible attachment site. Furthermore, straining in the low-flow velocity wedge region near the air-water-solid interface is a further retention mechanism for small sized particles when the AWI is not a favorable attachment site (Flury & Aramrak, 2017;Knappenberger et al., 2014).
X-ray microtomography (X-ray μCT) is a powerful tool to study phase distribution in porous media and to quantify the interfaces between the liquid transport medium, solid, and air. However, the spatial resolution of X-ray μCT is limited to some microns. This is why it is not possible to visualize the phase distribution within the pore space for fine textured soil. On the other hand, sand can be visualized at the pore scale and is an ideal porous medium to systematically study the transport processes under changing phase saturations, even though this implies a strong simplification compared to natural soil. μCT imaging of the phase distribution during flow experiments is a challenge as the phases and their interfaces have to be stable during the scan. To the best of our knowledge, the combination of X-ray μCT analysis of phase characteristics in porous media and unsaturated Ag-NP transport experiment has never been done so far.
Unsaturated column experiments are often controlled by a prescribed flux at the upper boundary and some negative pressure imposed at the lower boundary. This setup bears the risk to generate gradients in water content and water potential along the sample profile. Especially hydraulic state variables of coarse-grained materials are highly vulnerable to small changes at the boundaries, due to their highly nonlinear water retention characteristics and hydraulic conductivity. The Multi-Step-Flux (Kumahor, De Rooij, et al., 2015;Weller & Vogel, 2012) provides a setup to avoid these limitations. Our unsaturated experiments followed conditions studied for citrate-coated Ag-NPs (Kumahor, Hron, et al., 2015) and soil-aged Ag-NPs (Kumahor et al., 2016) for a better comparability.
In this study we used a cleaned quartz sand as model porous medium to study the impact of variable saturation on transport, retention, and remobilization of S-Ag-NPs. Specifically, we explore (i) the impact of presence or absence of NOM, (ii) the effect of different flow velocities, and, consequently, and (iii) the effect of different phase saturations and AWI. Two different transport models were used to describe the particle breakthrough and retention profile. The phase distributions (air, water, and solid) and the total interfacial surface (AWI and SWI) during water flow were analyzed via X-ray μCT.

S-Ag-NPs and NOM Preparation
A stock dispersion of citrate-coated silver nanoparticles (Cit-Ag-NPs) with Ag concentration of 1 mmol L −1 was synthesized applying citrate reduction method . The sulfidation of Cit-Ag-NPs was performed by rapid mixing 70.56 ml stock dispersion of Cit-Ag-NP with 1.44 ml stock solution of 100 mmol L −1 sodium hydrogen sulfide (NaHS, anhydrous, Alfa Aesar, ThermoFisher, Germany) in a 250 ml round-bottom flask under continuous stirring (stirring rate: 500 rpm) at room temperature. After 100 min, S-Ag-NPs were washed with a solution of 25 mg L −1 trisodium citrate (Sigma-Aldrich, USA) using centrifugal ultrafiltration membrane with a molecular weight cut-off of 50 kDa (Amicon Ultra-15 centrifugal filter device, Merck Millipore). Therewith the reaction was stopped, the residual, unreacted sulfide was removed, and NPs were stabilized. The total Ag concentration in prepared dispersion of S-Ag-NPs was 40.5 mg L −1 . The details on the preparation of S-Ag-NPs are given in Metreveli et al. (2020).
The intensity-weighted z average hydrodynamic diameter of S-Ag-NPs was measured using dynamic light scattering (DLS, laser wavelength: 658 nm, scattering angle: 165°, Delsa Nano C, Beckman Coulter, USA) after dilution (dilution factor: 1:3) of the stock dispersion of S-Ag-NPs in Milli-Q water (resistivity: 18.2 MΩ·cm, Direct-Q UV, Millipore, France). Size and occurrence of S-Ag-NPs were also characterized by transmission electron microscopy (TEM) and scanning TEM (STEM) with high-angle annular dark-field (HAADF) imaging (FEI Tecnai Osiris,accelerating voltage: 200 kV). In order to prepare the samples for TEM and STEM measurements, a drop of the stock dispersion was nebulized onto a copper grid (3 mm) covered with a combined holey and ultrathin (~3 nm) carbon film (Product No. 01824, Ted Pella, Inc., USA) applying an ultrasonic generator, which was developed at the Karlsruhe Institute of Technology. The distribution of silver and sulfide was determined by X-ray mapping using ChemiSTEM (FEI) system with four silicon drift detectors (SDD) of the Bruker Company.
A stock solution of Suwannee River NOM (SR-NOM, 2R101N, International Humic Substances Society) was prepared by dissolving 10 mg SR-NOM in 50 ml demineralized water. The solution pH was adjusted to 10 with 1 mol L −1 NaOH and filtered through a 0.2 μm cellulose nitrate membrane filter (Nalgene, USA). After stirring the filtrate with a magnetic stirrer for 24 hr, the stock solution was stored in a refrigerator.

Colloidal Stability of NPs and eDLVO Calculation
In order to evaluate the influence of the pH value and NOM on the particle size, aggregation dynamics, and ζ-potential of S-Ag-NPs, the NPs were dispersed in a solution of 1 mmol L −1 KNO 3 (ISO, Reag. Ph Eur, Merck, Germany) in the absence and presence of 5 mg L −1 SR-NOM and the pH values of 5, 6, and 7 were adjusted using HNO 3 (Suprapur, Merck, Germany) and NaOH (ROTIPURAN, p.a., Carl Roth, Germany). The samples with a total volume of 30 ml were prepared in plastic centrifuge tubes. The concentration of Ag in final dispersions was 1 mg L −1 . The samples with a pH value of 6 represented the S-Ag-NP dispersions applied in column experiments. The hydrodynamic particle diameter was measured with the same instrument as used for S-Ag-NPs diluted in Milli-Q water. For ζ-potential measurements, another instrument was applied (ZetasizerNano ZS, Malvern Instruments, UK). The determination of ζ-potential is based on the phase analysis light scattering (PALS) and mixed-mode measurement (M3) methods. The particle size and ζ-potentials were measured~5-10 min after sample preparation. The dynamics of S-Ag-NP aggregation was determined at several time points from 5 min to 4.1 days.
Using the measured values of NP size and ζ-potential the interactions between S-Ag-NPs with and without NOM aging and collectors representing SWI or AWI were calculated. For the calculations, the model of extended Derjaguin-Landau-Verwey-Overbeek (eDLVO) interactions was applied (Kumahor, Hron, et al., 2015). The detailed information regarding equations, physical constants, and input parameters used for the calculations is shown in the supporting information (SI).

Column Experiments
The transport experiments were conducted with a cleaned sand (Co. Sand-Schulz, Germany) packed to a bulk density of ρ b ¼ 1.52 g cm −3 . The grain size distribution was 0.1 to 0.3 mm with a median grain size of d 50 ¼ 0.2 mm (Kumahor, De Rooij, et al., 2015). A 1 mmol L −1 KNO 3 solution was used as background electrolyte solution (BES). Different transport and leaching solutions were prepared to test the impact of NOM: (i) dispersion of S-Ag-NP without NOM: 1 mg L −1 Ag in 1 mmol L −1 KNO 3 BES; (ii) dispersion of NOM aged S-Ag-NPs (NOM-S-Ag-NP): 1 mg L −1 Ag and 5 mg L −1 SR-NOM in 1 mmol L −1 KNO 3 BES; (iii) leaching solution without NPs and without NOM: 1 mmol L −1 KNO 3 BES; (iv) leaching solution without NPs and with NOM: 5 mg L −1 SR-NOM in 1 mmol L −1 KNO 3 BES. The mixtures were prepared 2 hr prior to each experiment to allow NP equilibration with SR-NOM and BES. Overall, eight transport experiments were performed with three different fluxes (j w * ¼ 17, 10, 2.5 cm hr −1 ) and saturations. Figure 1 gives an 10.1029/2020WR027074

Water Resources Research
overview about the different experimental setups, and Table 2 lists the information about the experimental conditions: the presence of NOM in transport and/or leaching solution, the mean pore water velocity, and the resulting water contents (θ [vol.%]) and pressure heads (h [cm]).
For each experiment a column was newly packed using the same procedure. There was only one replicate per experiment since the porous medium is well defined and macroscopically homogenous. This was confirmed by X-ray μCT measurements and the breakthrough curves (BTCs) of the tracer experiments. Further, the consistency of the results at different flow rates provided information on potential errors produced by the conditions not controlled in the experiment.
At first, a steady-state flow was established by applying 5 pore volumes (PV) BES at the predefined flux rate, followed by a tracer experiment with KBr (2.5 PV) and flushing with 3 PV of BES. The tracer experiment was conducted by adjusting the electrical conductivity (EC) of the tracer solution with KBr to 550 μS cm −1 and constantly measuring the EC of the outflow. The total application of transport solution containing S-Ag-NPs was then 8 PV, except for the 2.5 cm hr −1 experiment where 12 PV of transport solution was applied. The transport experiment ended after 8 PV flushing with BES (with and without 5 mg L −1 SR-NOM, respectively). At the end of the experiment, the sand column was dissected into 1 cm layers, dried at 105°C in an oven for 24 hr, and homogenized. The retained NPs were determined in a subsample of each layer. For the saturated experiments, the last PV was drained before slicing. During the experiments, pH and EC of the effluent were measured regularly. The samples of the outflows were stored in glass vials in a refrigerator at 4°C.
Saturated transport experiments ( Figure 1a) were conducted under two different flow conditions: 17 and 10 cm hr −1 . A cylinder made of polycarbonate (h ¼ 20 cm, d ¼ 9.4 cm) was placed on an airtight, funnel shaped tripod with a porous sintered glass plate (Degenkolb et al., 2019). To prevent air infusion at the lower boundary, a membrane (pore size ¼ 10 μm) was placed on top of the glass plate. The cylinder was filled with BES before adding the sand to a height of 10 cm. Entrapped air was removed by placing the setup inside a vacuum desiccator and degassing the water. For the experiment, both the setup with the soil column and the reservoir with the input solutions were placed on a balance to control the fluxes. A supernatant of 0.4 mm was adjusted to ensure water saturation throughout the experiment. A constant flux was applied at the lower boundary by a peristaltic pump, while a second pump provided the same flux at the sample top using an irrigation device with 21 equally distributed needles (inflow ¼ outflow). Gravimetric weights were automatically recorded every 60 s. The outflow was continuously sampled with a fractionator.
The unsaturated transport experiment ( Figure 1b) followed the protocol for Multi-Step-Transport experiments (Kumahor, De Rooij, et al., 2015), where the hydraulic conditions are controlled by the applied flux (Weller & Vogel, 2012). A tensiometer installed at 2.5 cm depth of the sample measures the water potential, which adapts to the imposed flux, while a pressure control system adjusts the pressure at the lower boundary to the water potential measured at 2.5 cm depth in real time. This procedure ensured unit gradient conditions along the sample profile, with constant soil matric potential and water flow driven by a unit gradient in gravity. The columns were placed on the same funnel shaped tripod with a porous sintered glass plate and membrane as described above. The dry sand was slowly filled in the column while tapping gently on the sides to ensure a uniform ρ b . The different transport and leaching solutions were applied with a peristaltic pump, and both the sample setup and reservoir were placed on a balance to control the fluxes and water contents gravimetrically. The outflow was continuously sampled with a fractionator.

Mathematical Modeling
For each experiment with known flux j w (cm hr −1 ), the basic transport parameters, mean velocity v (L T −1 ), and dispersivity D (L) were obtained by fitting the classical convection-dispersion model to the measured tracer breakthrough curves using CXTFIT (Toride et al., 1999). Model parameters for NP transport and retention were obtained by inverse modeling using HYDRUS 1D, V. 4.17 (Šimůnek & van Genuchten, 2008).
The applied models distinguish between two different attachment sites, which are supposed to interact differently with NPs. We used two different model approaches (i) the "Kinetic-Equilibrium" approach combines a kinetic attachment site representing SWI with an equilibrium attachment site described by a retardation factor representing AWI as proposed by Kumahor, Hron, et al. (2015) and (ii) the "Kinetic-Kinetic" approach assumes kinetic attachment at both interfaces using a modified version of the two kinetic sites model (Gargiulo et al., 2008;Schijven & Šimůnek, 2002). The model parameters were fitted to the experimental observations of both, the breakthrough curve and the retention profile simultaneously.
Kinetic-Equilibrium model. The classical convection-dispersion equation is extended by a retardation factor R (-) describing reversible interaction while the attached mass of NPs is in equilibrium with the mass of NPs remained in aqueous phase. Additionally, a second type of attachment site is described by a nonequilibrium kinetic process described by an attachment, K att (T −1 ), and detachment rate, K det (T −1 ). The transport of NPs is then described as where t (T) is time and z (L) is the spatial coordinate. The dimensionless retardation factor is given by Assuming linear sorption (Abraham, 2015), R relates to an equilibrium attachment coefficient, K eq : where ρ b (M L −3 ) is the soil bulk density, θ is the volumetric water content (-), S NP (M M −1 ) is the mass of attached NPs relative to the mass of solid phase, and C (M L −3 ) is the concentration of NPs in the aqueous phase. The nonequilibrium reaction term is given by The time-and depth-dependent dimensionless straining function ψ (-) is defined as (Bradford et al., 2006) where S NP, max (M M −1 ) is maximum mass of attached NPs relative to the mass of solid phase, d 50 (L) is the median diameter of sand grains, and β (-) is an empirical parameter, which describes the shape of the retention profile along the flow path. An equivalent model approach was used by Torkzaban et al. (2008) to describe the transport of colloids in unsaturated sand.
Kinetic-Kinetic model. This model assigns the retention of NPs to two different reaction sites (S 1 ¼ reversible and S 2 ¼ irreversible) and hence divides the retention into two different kinetic processes (Hassanizadeh & Schijven, 2000). The transport equation is defined as In analogy to Equation 4, the first kinetic site (S 1 ) provides the possibility of NP attachment and detachment (reversible): with a depth-dependent straining function ψ S1 (-) (Bradford et al., 2003) Attachment to the second site (S 2 ) is irreversible and defined as Here, ψ S2 is a time dependent retention function given as The Kinetic-Kinetic model of the BTC and the measured data of the retention profiles were used to calculate the mass balances.

Imaging With X-Ray μCT
The transport experiments provided the pressure heads and water contents for the different applied flow rates. To simulate the transport condition on a smaller scale, a flow tube (h ¼ 8 cm, d ¼ 1.5 cm) made of glass was packed with sand having the same bulk density and placed inside an X-ray μCT device (XTH225, Nikon Metrology, Belgium). Water flux was provided by a single dripper to the sample top. At the lower boundary a hanging water column was connected to a sintered glass plate ( Figure 1c). Water pressure was adjusted by the hanging water column to the same value as measures for the same flow rate during the transport experiments in the larger columns (Table 2). Before μCT scanning, the flux was equilibrated for 4 PV. Subsequently, to reduce possible pulsation of AWI, the connection to the atmosphere was closed and water flow was equilibrated for two more PVs. During constant irrigation at a defined flow rate and water potential, the setup was scanned for 16 min at a resolution of 7.5 μm per voxel with 1,200 projections. Starting from capillary saturation the flow rate and potential at the lower boundary were stepwise reduced to 17, 10, and 2.5 cm hr −1 .
Image segmentation into the three phases (solid, air, and water) was done using the software packages FIJI ImageJ (Schindelin et al., 2012) and QuantIm (Vogel et al., 2010). The reconstructed images were filtered with a 2-D-nonlocal means filter to reduce image noise (Buades et al., 2005). Thereafter, the gray values 10.1029/2020WR027074

Water Resources Research
were segmented into the three phases using the minimum error thresholding method. Partial volume effects were reduced by a majority filter (Schlüter et al., 2014). The phases distributions of solid, air, and water and their interfaces were determined for three subvolumes (750 * 750 * 370 voxel ¼ 0.086 cm 3 ) with the MorphoLibJ plugin for FIJI ImageJ (Legland et al., 2016) to capture potential spatial heterogeneities. The interfacial area A ab between two phases a and b within the three-phase system (a, b, and c) was calculated based on the total surfaces of the single phases A x as (Schlüter et al., 2014;Vogel et al., 2010) The determined volumetric phase contents and their interfacial areas were upscaled to the larger 10 cm columns (694 cm 3 ). The ratio of Geodesic distance and Euclidean distance (cm cm −1 ) from the image bottom to the top was calculated for the segmented water phase as a measure for its tortuosity.

Ag Quantification
Concentrations of Ag-containing NPs in dispersions were measured by single particle ICP-sector field mass spectrometry (Element XR, Thermo Fisher) with 3 ms integration time and analyzed with a modified version of SPC tool developed by RIKILT-Wageningen UR. The Ag standard curve (0.02 to 2 ng ml −1 ) was calibrated on ICP-Multi-element solution VI (Merck CertiPur, Germany) and ICP-MS-Multi-element solution 2 (SPEX CertiPrep, USA). For quality control, a Pt-nanoparticle standard (50 nm, citrate coated, NanoComposix, Czech Republic) was used. In order to determine total Ag concentration in aqueous samples, 5 ml HNO 3 (65%, suprapur, Merck, Germany) was mixed with 5 ml sample followed by digestion in a microwave system (200°C, 30 min). Digested samples were diluted with Milli-Q-water according to the expected concentrations for the final measurement. S-Ag-NPs retained in sand were digested by adding 5 ml HNO 3 (65%, suprapur, Merck, Germany) and 0.5 ml H 2 O 2 (30%, suprapur, Merck, Germany) to 0.6 g sand (dried at 105°C for 24 hr, and homogenized) followed by microwave digestion (200°C, 30 min). After separation of residual sand and diluting with Milli-Q-water, the Ag concentration in supernatant was measured by ICP-Quadrupole mass spectrometry (IcapQs, Thermo Fisher, Germany). Finally, the Ag content in sand was calculated from measured Ag concentration and sand mass.

S-Ag-NP Characteristics
In the absence of NOM, the mean hydrodynamic diameter of S-Ag-NPs dispersed in 1 mmol L −1 KNO 3 solution at pH 6 corresponding to the chemical conditions in transport experiments was 60.5 ± 1.4 nm ( Figure S1a, SI). A higher pH value of 7 did not lead to the significant changes in particle size. Only at lower pH of 5, S-Ag-NPs aggregated to a hydrodynamic diameter of 95.9 ± 4.1 nm~10 min after sample preparation in the absence of NOM. In the presence of SR-NOM, the hydrodynamic diameter of S-Ag-NPs was in the range of 58-60 nm at all pH values indicating high colloidal stabilization of the NPs by SR-NOM. Similar tendencies were observed for particle size distribution ( Figure S1b, SI). Interestingly, the size distribution of S-Ag-NPs determined by DLS in the absence as well as in presence of SR-NOM was broad, covering hydrodynamic diameters from about 15 nm to 500 nm. The ζ-potential determined for S-Ag-NPs at pH 6 and 7 in the absence of NOM was negative in the range from −33 to −36 mV ( Figure S1c, SI). The reduction of pH value to 5 resulted in reduced negative value of the ζ-potential (−25 mV). In the presence of SR-NOM the ζpotential of S-Ag-NPs was more negative than in the absence of NOM. Furthermore, in the presence of SR-NOM, the ζ-potential decreased with increasing pH value from −36 to −64 mV. At pH values between 6 and 7, the hydrodynamic diameter of the NPs in the absence as well as in the presence of SR-NOM did not change significantly within the time period of 4.1 days and remained in the range of initial nanoparticle size ( Figure S1d, SI). The TEM measurements showed that the S-Ag-NPs in stock dispersion were present as single sulfidized Ag-NPs and their small aggregates (Figures 2a and 2b). Furthermore, the NPs were partly or completely sulfidized (Figures 2c and 2d, and S2a and S2b, SI).
The dimensionless energy profiles calculated according to the eDLVO theory for the interactions between unaged or NOM-aged S-Ag-NPs and collectors representing SWI or AWI are shown in Figure S3 (SI). Calculations were performed for the pH values of 5, 6, and 7. Since the information on the contact angle 10.1029/2020WR027074

Water Resources Research
is not available for S-Ag-NPs, we used a contact angle of 40°reported by Hu et al. (2004) for citrate-coated Ag-NPs as well as slightly lower (30°) and slightly higher (50°) values. Depending on the pH values and the presence or absence of NOM, the calculations showed an absence of secondary minimum and high energy barriers ranging from 75 to 144 kT at the separation distances of 1.0-1.6 nm for NP-SWI interactions. The changing of hydrophobic properties by variation of contact angle of NPs did not lead to changes in NP-SWI interaction energy profiles. This is in contrast to NP-AWI interactions. Here, the high of energy barriers ranged from 3.4 kT (S-Ag-NPs, pH 5, contact angle: 40°) to 121 kT (NOM-S-Ag-NPs, pH 7, contact angle: 30°) at the separation distances of 2.6-9.0 nm and decreased or completely disappeared with increased contact angle. For the contact angle of 50°, the energy barrier with a high of 3.6 and 29 kT was only observed in the case of NOM-S-Ag-NPs at pH values of 6 and 7, respectively. In contrast to the NP-SWI interactions, a secondary minimum was observed in most cases for the NP-AWI interactions. The depth of the secondary minimum ranged from −2.1 to −11.2 kT at the separation distances of 35-62 nm.

Transport Experiments
The Multi-Step-Flux setup enabled the study of S-Ag-NP and NOM-S-Ag-NP transport under well-defined conditions. The achieved water contents, pressure heads and the experimentally realized transport velocities for the different experiments are shown in Table 2. The low standard deviations of the hydraulic properties document the stable experimental conditions. A reduction in fluxes (j w * ¼ 17, 10, and 2.5 cm hr −1 ) resulted in reduced water contents (θ ¼ 30, 27, and 19 vol.%) and decreasing pressure heads ( Table 2). The determined mass balance for each experiment is listed in Table 1, resulting in a total recovery rate ranging from 80.7% to 106.3%.  experiments with Cit-Ag-NPs at pH 5 for identical conditions in terms of sand, hydraulic conditions, and input concentration of Ag and BES as conducted by Kumahor, Hron, et al. (2015). Table 2 provides all fitted parameters and their standard error coefficients for the two different models. For the Kinetic-Equilibrium model the variables S max , k att , k det , k eq , and β were fitted, for the Kinetic-Kinetic model the variables S max , k att1 , k det1 , k att2 , and β.
Overall, the transport of S-Ag-NPs through the sand columns was changing with the presence of NOM, the applied flux, and the water content. In the absence of NOM in transport solution (unsaturated, high flow scenario in Figure 3a, gray dots) no particle breakthrough was observed and nearly 65% of the applied NPs were attached to the upper 2 cm of the column. By contrast, in the presence of NOM, particle breakthrough was observed for all experiments but varying in time and total amount. Comparable results in BTCs and

10.1029/2020WR027074
Water Resources Research retention profiles were found for the saturated experiments and the unsaturated high flux scenarios. They were marked by a moderate increase of particle breakthrough up to 65% of the initial concentration and a sharp decrease after starting the leaching. The four high and medium experiments differed in their starting points of particle breakthrough. The first measured particle breakthrough (Ag ≥ 0.01 mg L −1 ) occurred around 2.5 PV for the high flow scenarios with NOM-S-Ag-NP and after 3.4 PV for the saturated intermediate flow scenario, respectively. Furthermore, a tailing of Ag output was observed for both unsaturated experiments; that is, concentrations did not drop back to 0. For the unsaturated experiments, a decrease in water flux from 10 to 2.5 cm hr −1 and concurrent reduction in water content from 27 to 20 vol.% caused a strong decrease of the total particle breakthrough, as revealed by the lower peak of the BTC, while the particle retention was increased. The initial particle breakthrough occurred after 4.3 and 4.7 PV at 10 cm hr −1 and after 5.2 PV at 2.5 cm hr −1 . All NOM-S-Ag-NP experiments showed a sharp decrease in particle breakthrough shortly after changing to the leaching solutions. The presence or absence of NOM in the leaching solution had no pronounced effect on the BTCs and profiles.
The NOM-S-Ag-NP experiments received 7 PVs for high and intermediate flows, and 12 PV for the low-flow scenario, before switching back to the leaching solution. The breakthrough at 7 PV was compared between the scenarios, since they all received constant application of NOM-S-Ag-NP up to this point. The percentage of particle breakthrough at 7 PV was highest for the saturated (17 cm hr −1 : 13.8% and 10 cm hr −1 : 9.52%) and unsaturated high flow experiments (11.3% and 12.4%), whereas breakthrough was reduced for the intermediate (4.3% and 3.9%) and more so for the low-flow scenario (0.8%). Consequently, the percentage of total Ag breakthrough (Table 1) was highest for the saturated (32.3%) and unsaturated (24.5% and 28.6%) high flow scenario and the saturated intermediate flux (24.0%). For unsaturated conditions, the particle breakthrough decreased to 6.3% and 8.2% at 10 cm hr −1 , and to 6.4% at 2.5 cm hr −1 .
The retention profiles of NOM-S-Ag-NP experiments at 17 and 10 cm hr −1 followed the same general pattern of a strong NP attachment (S/C 0 ratios of 3 to 8) close to the column surface, a sharp decrease, and a small but rather constant attachment over the sample depth from 3 to 10 cm. While the majority of these retention profiles can be defined as hyperexponential, some of them are nonmonotonic (Goldberg et al., 2014). The S-Ag-NP experiment (gray dots) showed the same behavior, but in a more pronounced form leading to high retention at the sample top and reduced concentrations below 6 cm depth. For the 2.5 cm hr −1 experiment, the retention profile was non-monotonic, where the concentration peaked in the first 3 cm (S/C 0 ratios of 1.9 to 3.6) before decreasing with depth.
In comparison with the BTC of Cit-Ag-NPs, the unsaturated high and intermediate flow experiments with NOM-S-Ag-NPs showed an initial particle breakthrough that was delayed by around 2 PV. The retention profiles, on the other hand, had similar general shapes, but the concentrations in the lower part of the samples were systematically higher for NOM-S-Ag-NPs. For the low-flow experiments Cit-Ag-NPs did not show any breakthrough in contrast to NOM-S-Ag-NPs.
From visual inspection, the Kinetic-Kinetic model resulted in the best fits of the tested models to describe both types of information, the BTCs and retention profiles (Figure 3 solid lines and Table 2, b). Two kinetic processes were superior to describe the slope of particle breakthrough of the high and intermediate flow conditions, as well as the sharp decrease in NP retention below the column top. For the low-flow scenario, no differences between the model qualities to describe the BTC were observed. Overall, the Kinetic-Kinetic model resulted in a good approximation of the data as also suggested by R 2 ¼ 0.79 to 1.00. In some cases, the R 2 of the Kinetic-Equilibrium model indicated a better fit, which was due to the high impact of the uppermost concentration value in the retention profile on the calculated R 2 value. Here, the hyperexponential form of the Kinetic-Kinetic model made the fit highly sensitive to the sharp drop between the top value and the spatially continuous values below, especially if the model was not able to fit both parts accurately. These differences then caused a higher root-mean-square-error for the Kinetic-Kinetic model, even though the rest of the retention profile and the BTC was fitted better. It is noted that the measured concentrations of the retention profile were averaged over 1 cm depth intervals of mixed material, which may lead to deviations between model results and measurements especially for the highly nonlinear profiles.
Model fitting resulted in large standard error coefficients and high correlations between parameters. Hence, a physical interpretation of the fitted parameters is hardly possible and the applied models are merely a description of the measurements. Nevertheless, the Kinetic-Kinetic model allows for separating the 10.1029/2020WR027074

Water Resources Research
retention profiles into two components that can be related to two different attachment sites. While Site 1 shaped the hyperexponential rise near the top, Site 2 created the low but rather uniform attachment over the sample depth ( Figure S3, SI). The two different sites have been assigned to a SWI and an AWI, respectively (Kumahor et al., 2016). The best model results were obtained when attachment parameters were fitted for two sites, even though for the saturated samples only erratic occurrence of an air phase due to entrapment would have been expected. As a result, and due to parameter correlations, there was no obvious trend for single parameters of both models to be related to the applied fluxes, the changing transport velocities and water contents, and hence to be related to the importance of the AWI.

Characterization of the Phase Distribution in Porous Media
The image analysis of the flow cell for the three different fluxes resulted in approximately the same water contents as measured in the larger column experiments (Tables 3 and 2, respectively). Based on the μCT images it was possible to visualize the phase distributions of solid, water, and air, and to quantify their interfaces. Figure 4 shows the phase distributions for a representative volume (150 3 voxel) for all applied fluxes. The depth profiles of the single-phase distribution confirmed the establishment of a steady-state flow ( Figure S4, SI). It is shown that the analyzed part of the flow tube (Figure 4, left, and Figure S4, SI) was Note. Phases distributions of water, air, and solid, their interfaces (upscaled to the total interfacial area of the 10 cm column), and the tortuosity of the water phase are presented as mean and their standard deviation of three subvolumes. The tortuosity was determined as the ratio of Geodesic distance (cm) to the Euclidean distance (cm).

Water Resources Research
homogenously packed, and that the internal structure was macroscopically homogenous, a prerequisite for upscaling and interpretation of the determined data.
The total interfacial areas summarized in Table 3 were upscaled to the volume of the large column (697 cm 3 ). The AWI for the saturated sample can be attributed to air entrapment. In contrast to the larger sample, the flow cell was not degassed before the start of the experiment. Overall, the differences in SWI between the scenarios were less pronounced than the differences in AWI. While the AWI was increasing from 0.9 to 6.2 m 2 , the SWI was decreasing from 12.6 m 2 to 11.0 m 2 . The reduction in SWI was an artifact of image resolution, as a thin water film should have remained at the solid surfaces. Hence, the SWI was still the dominant interface. The main jump in AWI increase was determined for the flux step from 17 cm hr −1 and 31 vol.% water content to 10.0 cm hr −1 and 27 vol.% water content. Here, the area of AWI was increasing by a factor of 2.5. Further, with decreasing water content the tortuosity was increasing. The mean transport distance along the path way was about 6.4% longer for the saturated experiments and about 12.9% for the unsaturated low-flow scenario. In comparison to the increase in AWI, the impact of the decreasing water content on tortuosity was less pronounced.
The measured phase distributions and interfacial areas were strongly correlated to the measured characteristics of the BTCs (Figures 5 and S6, SI). The increase in AWI was related to the decreasing particle breakthrough with decreasing flux and, hence, decreasing water content. Figure 5 shows the three described characteristic parts of the NOM-S-Ag-Np BTCs as a function of mean transport velocity (a) and the AWI (b): the delayed initial particle breakthrough (left), the decreasing amount of particle breakthrough after 7 PV of NP application (center), and the total reduction in particle breakthrough after 16 PV (right). The blue arrows mark the shifts in particle breakthrough for the experiments with the same fluxes (j w * ¼ 17 and Figure 5. The determined initial particle breakthrough (PV), the cumulative NOM-S-Ag-NP particle breakthrough after 7 PV (percentage of particle mass application), and the total particle breakthrough (percentage of particle mass application). Illustrated as a function of (a) the mean transport velocity and (b) the determined area of air-water interfaces (AWI). The blue arrows mark the shifts for experiments with comparable fluxes, j w , but different saturations.
10 cm hr −1 ) but different saturations (black ¼ saturated, gray ¼ unsaturated). The total breakthrough of the 2.5 cm hr −1 experiment is not included, due to the additional application of 4 PV particles.
The regressions suggest that the initial particle breakthrough (PV) of the different experiments was affected by both, the AWI and the transport velocity. The decrease in transport velocity could explain 65% of their variance, the increasing AWI 86%. For the 17 cm hr −1 experiments (saturated and unsaturated) with small differences in their mean transport velocity and AWIs, the initial particle breakthrough occurred around the same PV, while the differences between the 10 cm hr −1 experiments became much larger ( Figure 5, left). In fact, the differences between the saturated samples with different fluxes were even smaller than the differences between experiments with the same fluxes but different saturations, indicating that the initial BT was largely affected by the additional air phase and not the decreasing velocity. Compared to all other measured parameters (tortuosity, saturation, SWI, Figure S7, SI), it is shown that with consecutive particle breakthrough the total interfacial area provided the best explanatory variable, which was mainly dependent on the increase in AWI.
It is noted that for the unsaturated experiments, the different soil hydraulic characteristics are interrelated: a reduction of applied fluxes caused reduced water contents and hence increased AWIs and tortuosity of the water pathway. Further, the determined transport characteristics are strongly correlated, as a delayed particle breakthrough resulted systematically in a reduced total particle breakthrough.

Discussion
The transformation of Cit-Ag-NPs to S-Ag-NPs and NOM-S-Ag-NPs changed their physical and chemical properties. The changes in morphology, structure, and composition as a result of sulfidation process was also observed in Metreveli et al. (2020) showing that the S-Ag-NP aggregates were composed of few primary Ag-NPs sulfidized and "fused" together by silver sulfide (Ag 2 S) nanobridges. In addition, Ag 0 -Ag 2 S core-shell structures, asymmetrical sulfidation, and hollow regions for S-Ag-NPs were observed. The changes in surface properties of Ag-NPs during sulfidation may have an influence on their mobility in porous media. The surface of the silver sulfide, which tends to be originally hydrophobic can interact with the hydrophobic chain of organic surfactants (Fuerstenau, 1980). Consequently, despite the small amount of citrate bound onto the surface of NPs, an increased interaction with AWI which is known as super hydrophobic (Hu et al., 2004) is expected for S-Ag-NPs. Furthermore, the DLS measurements showed that the sulfidation resulted in increasing hydrodynamic diameter of NPs from~30 to 60-80 nm. The SEM and TEM measurements revealed an increasing size during sulfidation, from 48.8 ± 20.1 nm for Cit-Ag-NPs to 81.1 ± 52.2 nm S-Ag-NPs (Metreveli et al., 2020).
The aggregation of S-Ag-NPs in the absence of NOM at pH 5 can be explained by reduced negative value of the ζ-potential and thus by reduced electrostatic repulsion forces between NPs as a result of increased protonation of their surface. The aging of S-Ag-NPs in the presence of NOM most probably leads to the NOM coating on the NP surface suggested by increased negative ζ-potential. As a consequence, S-Ag-NPs became more colloidally stable due to the increasing electrostatic repulsion forces between NPs. The high colloidal stability of S-Ag-NPs at pH values of 6-7 observed in the absence and presence of SR-NOM within the time period of~4 days suggests that the NPs did not aggregate during the transport experiments.
The BTCs and retention profiles of NOM-aged and non-aged S-Ag-NPs in quartz sand differed to previous experiments with stabilized Ag-NPs. The main difference was the shift in time for initial NOM-S-Ag-NPs breakthrough of 2.5 to 5.2 PV. Most studies reported for stabilized Ag-NPs an earlier breakthrough at around 1 PV independent of the saturation, for example, for saturated experiments with and without NOM (El Badawy et al., 2013;Kanel et al., 2015;Taghavy et al., 2013), and unsaturated experiments (Liang, Bradford, Šimůnek, Heggen, et al., 2013;Yecheskel et al., 2016). Under the same unsaturated hydraulic conditions, Kumahor et al. (2016) observed for Cit-Ag-NPs and Cit-Ag-NPs aged with soil solution an initial breakthrough between 0.5 and 1.7 PV. Only for the low-flow scenario with Cit-Ag-NPs a comparable retardation was observed at a solution pH of 9, while at pH 5 no particle breakthrough occurred. The retardation of particles was assigned to a reversible attachment process at the AWI, which is not transferable to our data as the saturated experiments also showed particle retardation. A reduction of total Ag-NPs breakthrough and a reduced release rate have been observed in previous studies when potential 10.1029/2020WR027074

Water Resources Research
collector sites were present (El Badawy et al., 2013;Taghavy et al., 2013) and which was further reduced with decreasing transport velocity (Hou et al., 2017). Both, the retardation and the reduction of breakthrough for the saturated experiments with NOM-S-Ag-NPs suggest an attachment process at the solid phase. Consequently, the factor transport velocity cannot explain the further increase in particle retention of the unsaturated experiments with decreasing fluxes. On the other hand, the energy barriers determined by eDLVO calculations for NP-SWI interactions were much higher (75-144 kT) than the energies of Brownian diffusion (20 kT) (Elimelech et al., 1995) significantly reducing the probability for the irreversible attachment to SWI in primary energy minimum. The eDLVO calculations do not consider the complex shapes and morphology of S-Ag-NPs (Metreveli et al., 2020) as well as the nanoparticle attachment to SWI due to the complex formation between NOM and possible cations present on the sand surface. Furthermore, the presence of small amount of AWI in saturated column despite the degassing of water prior to the experiment cannot be fully excluded providing additional attachment places for NPs. At saturation, the lower transport velocity only slightly affected the BTC.
The complete particle retention of S-Ag-NP without NOM aging in our study was in contrast to the findings of Yecheskel et al. (2016) where sulfidized Cit-Ag-NPs behaved like a conservative tracer in sand. The authors used a different sulfidation method to generate the NPs where a higher concentration of citrate was used, the sulfidation was performed by granular pyrite, and no NP washing of surplus citrate and sulfur was described.
The sulfidised NPs only slightly reduced their ζ-potential (−55.5 to −53.3 mV), and their shape was round compared to the elongated S-Ag-NPs used in this study. Differences in S-Ag-NP production have also been used in the study by Zhu et al. (2016) who investigated the effect of particle aging on flocculation and filtration of PVP stabilized Ag NPs and found a reduced stability in self-aggregation and increased mobility for sulfidized NPs. The high of the energy barrier obtained by eDLVO calculations at pH 6, reflecting on average the pH range in the input solutions (pH 6.3 ± 0.2) and effluents (pH 5.6 ± 0.3) of the transport experiments, was 23.3 kT for the interactions between AWI and unaged S-Ag-NPs with a contact angle of 40°and decreased with decreasing pH value and completely disappeared with increasing contact angle. Since the energies induced by Brownian diffusion motion are in the range of 20 kT (Elimelech et al., 1995), the unaged S-Ag-NPs can easily overcome the energy barriers <20 kT and irreversibly attached to AWI at the primary energy minimum. At pH values of 5.3-5.9 measured in the effluents, the probability for the irreversible attachment of unaged S-Ag-NPs to AWI is significantly increased since energy barriers smaller than 20 kT are expected. Since for silver sulfide surface generally high hydrophobicity is suggested (Fuerstenau, 1980), we cannot exclude that unaged S-Ag-NPs have a higher contact angle of, for example, 50°at which the energy maximum is completely disappeared and the irreversible attachment to AWI is expected. While it became obvious that additional research is required to study the role of different sulfidized forms of NPs and their surface hydrophobicity on the fate and transport of chemically transformed and aged Ag-NPs in porous media, our study suggests that NOM aging may reestablish mobility of otherwise immobile NPs.
Hyperexponential and nonmonotonic retention profiles have been observed for Ag-NPs at various saturation, flow rates, and solution conditions (Degenkolb et al., 2019;Kumahor et al., 2016;Liang, Bradford, Šimůnek, Heggen, et al., 2013). An excessive deposition of NPs at the column inlet have been assigned to chemical particle heterogeneity or NP aggregation (Goldberg et al., 2014;Yecheskel et al., 2016), to different hydrodynamics (Bradford et al., 2011) and straining in the critical zone of the uppermost mm (Kumahor, Hron, et al., 2015). In our study, only the S-Ag-NP experiment provided an ideal hyperexponential retention profile. The profiles of NOM-S-Ag-NPs became slightly nonmonotonic. The saturated and the high flow experiment with NOM in the leaching solution, hence the most favorable conditions for NP mobility, resulted in more uniform transport into depth, which is in accordance with Yecheskel et al. (2018). The application of NP solution by droplets produced locally an increased concentration, which might support attachment to the AWI at the surface of an unsaturated sand column. Nevertheless, also for the saturated samples an increase of particle retention at the uppermost sample was observed even though a continuous ponding water table of 4 mm should avoid the establishment of a critical zone.
The Kinetic-Equilibrium model, which was used successfully in previous studies (Kumahor et al., 2016;Kumahor, Hron, et al., 2015), was able to fit the BTC and retention profiles only for some of the different scenarios. The use of a two Kinetic sites model with depth-dependent attachment resulted in a reasonable description of all experimental data, except for nonmonotonous features, which are not captured by any existing model so far. Furthermore, the usage of a kinetically controlled NP attachment instead of an 10.1029/2020WR027074 Water Resources Research equilibrium based process was highly recommended by Praetorius et al. (2014) as equilibrium principles are not applicable to NPs due to their thermodynamically unstable occurrence in the environment. In contrast to previous studies where the shape factor β was fixed for the tested porous media (sand β ¼ 0.43, e.g., Bradford et al. (2003)), β was a free parameter to give the model more flexibility in describing the different retention profiles (Goldberg et al., 2014). Unfortunately, a quantitative physically based characterization of two different attachment sites, that is, the increasing impact of an air phase, is not available. In addition, with the high correlations between the fitting parameters and their standard errors, a clear assignment of the interfaces to the two kinetic model sites was not possible.
The detection of phase characteristics in the porous media during flow conditions via X-ray μCT enabled a detailed analysis of the different drivers of particle retention. The observed BTC characteristics were related to the interfacial areas of the water phase to solid and air, which confirms the model results that two collectors contributed to the retention of NOM-S-Ag-NPs. It was shown that both the tortuosity of the mobile phase and the SWI were only slightly affected by the decreasing flow rate. Even though the SWI provided quantitatively the dominant interface and was involved in the particle attachment, the determined changes in BTC characteristics were highly correlated to the disproportionate increase of AWI with decreasing flux and water content. Additionally, a larger velocity contributes to earlier initial breakthrough and higher mass loss of particles in effluent after 7 PV.
For colloid transport in unsaturated media it was shown that the flow rate affected the retention of colloids retained in the secondary energy minimum but that the reduction in water content dominated the reduction of particle breakthrough (Knappenberger et al., 2014). Considering only electrostatic attraction forces, Wan and Tokunaga (2002) demonstrated in bubble column experiments that the negatively charged AWI interacts only with positively charged particles. However, for the SWI the classical DLVO calculations predicted that due to a secondary energy minimum NPs can still attach to quartz by van der Waals forces . The interplay between electrostatic and van der Waals forces is dependent on the ionic strength of the aqueous phase, electrical charge of particles and collectors, and their surface composition. Furthermore, it was also shown that hydrophobicity, which is not covered by classical DLVO theory can highly influence the interactions between particles and collector and that the AWI is the favorable attachment site when hydrophobic interactions are considered (Kumahor et al., 2016). According to the eDLVO calculations performed in our study, the irreversible attachment to AWI in primary energy minimum is predominantly expected for unaged S-Ag-NPs and cannot be excluded for NOM-aged S-Ag-NPs especially at low pH values and high contact angles of NPs. Due to the higher energy barriers as a result of more negative ζ-potentials compared to unaged S-Ag-NPs, predominantly the reversible attachment to AWI in secondary energy minimum can be suggested for NOM-aged S-Ag-NPs. The calculations using eDLVO theory in the earlier study showed that in contrast to the SWI, the citrate coated as well as soil-aged Ag-NPs can also attach to the AWI in secondary energy minimum (Kumahor et al., 2016;Kumahor, Hron, et al., 2015). The aging of particles with NOM can change the hydrophobicity and wettability of particle surfaces. High hydrophobicity of NOM-S-Ag-NPs can lead to the increasing probability for the attachment of NPs on AWI in the secondary energy minimum. Since attachment in the secondary energy minimum is often reversible, the NPs attached on AWI can easily be remobilized. On the other hand in addition to the possible irreversible attachment to AWI in primary energy minimum, the increased retention of NPs with decreased water content can be explained by straining on the wedge region between air-water-solid interface due to capillary pinning and the local low-flow velocity (Flury & Aramrak, 2017). A better link between CT-derived phase characteristics and their translation into a physically based transport model description would help gaining new insights into the processes governing NP retention.

Conclusion
We could demonstrate that the transformation of engineered Ag-NPs to S-Ag-NPs and, further, to NOM-S-Ag-NPs, which is the main expected transformation process considering real emission pathways and natural environmental conditions, has an impact on NP mobility in porous media. While the sulfidation of the particles induced a complete retention in quartz sand, the aging with NOM reestablished a certain mobility. As particles did not change their size nor their aggregation behavior due to the presence of NOM, it seemed that the physicochemical transformation of the particle surface caused the increase in mobility.
In this study, we used a combination of transport experiments at different saturations, numerical models, and X-ray μCT characterization of the porous media to systematically study the impact of interfacial areas and flow velocity on NOM-S-Ag-NP mobility in porous media. Based on the results presented herein, it was shown that 1. X-ray μCT imaging of the phase distributions during stationary flow was a powerful setup to quantify the extent of potential attachment sites in unsaturated porous media. There was a strong link between transport behavior of NPs and the structure of the phase distribution including the different interfaces. 2. A two-kinetic sites model with different retention characteristics was a suitable model to represent the sensitivity of the transport behavior to the hydraulic conditions as a function of the flow boundary conditions. 3. The retardation of particles was mainly dependent on the interfacial areas between the solid-water and air-water phase. The decreasing amount of particle breakthrough at reduced saturation was dependent on both, the increasing AWI and the decreasing flow velocity. 4. The transformation of Cit-Ag-NPs to NOM-S-Ag-NPs, as it is typically expected when released into the environment, decreases the mobility of these particles.

Data Availability Statement
All data are archived and provided online by the UFZ: https://www.ufz.de/record/dmp/archive/8411/en/.