The effect of flow on radiolysis in liquid phase-TEM flow cells

Applying a continuous flow to rinse radiolytic species from the irradiated volume is a widely proposed strategy to reduce beam-related artefacts in Liquid-Phase Transmission Electron Microscopy (LP-TEM). However, this has not been verified either experimentally or theoretically to date. Here we explore an extended numerical model implementing radiolytic chemistry, diffusion and liquid convection to study the peculiarities of beam-induced chemistry in the presence of a flowing liquid within a heterogenously irradiated nanoconfined channel corresponding to a LP-TEM flow cell. Intruigingly, the concentration of some principal chemical species, predominantly hydrogen radicals and hydrated electrons, is found to grow significantly rather than to decrease in respect to zero-flow when moderate flow conditions are applied. This counterintuitive behaviour is discussed in terms of reactants’ lifetimes, spatial separation of the reaction network and self-scavenging by secondary radiolytic species. In the presence of a flow the consumption of highly reactive species is suppressed due to removal of the self-scavengers, and as a result their concentration in the irradiated area increases. A proof of concept for the supply of scavengers by the flow is demonstrated. Unravelling the effect of flow on radiolysis spawns direct implications for LP-TEM flow experiments providing yet one more control parameter for adjusting the chemistry in the irradiated/imaging area, in particular for mitigation strategies by continuous supply of scavengers.


Introduction
Liquid-Phase TEM and radiolysis Liquid-Phase Transmission Electron Microscopy (LP-TEM) has evolved into a convenient experimental technique for the imaging of dynamic processes (biological, chemical, electro-chemical, etc.) at the nano-meter scale in wet-chemical environments [1][2][3]. Besides the uncontested benefits, such as high spatial resolution of TEM [4], utilising electrons for imaging caused several challenges to arise [5][6][7][8]. The main issue of using highenergy electrons for imaging is the inevitable creation of radiation damage in the irradiated matter, with radiolysis of water molecules being the most relevant in LP-TEM [3].
Briefly, radiolysis describes the scission of molecules (most importantly water) upon high-energy ionizing irradiation resulting in the creation of molecular, ionic and radical species [9][10][11][12]. The latter ones are highly reactive, thus can interfere with the (bio-) chemical processes of interest [6,13].
Since the early days of LP-TEM, radiolytic effects have been ambiguously discussed in the community: when the technique was limited to static liquid cells (LCs; e --beam-transparent reservoirs entirely enclosing the sample solution, thus protecting it from the surrounding vacuum of the microscope), the e --beam itself was the only (thus well-appreciated) stimulus to induce nanoscale dynamics [14][15][16]. With the evolution of LCs, a growing Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. number of stimuli (e.g., temperature, electric bias, solution replacement and mixing) were incorporated thus enabling the reproduction of macroscopic (electro-)chemical experiments at the nanoscale within the TEM [17][18][19][20]. Consequently, radiolytic effects have been increasingly perceived as obstructive artefacts that require thorough consideration and mitigation [6,20].

Strategies to account for radiolysis
The strategies developed to deal with radiolytic effects are twofold: (1) a profound theoretical understanding & description of the radiolytic processes in complex reaction media helping to explain experimental observations and (2) the mitigation of the radiolytic species in the irradiated volume and limiting their effects on physicochemical processes under imaging [21].
The field of numeric modelling for theoretical description of radiolysis in LP-TEM was initiated predominantly by the work of Schneider and co-workers [13] who adopted existing knowledge [12] in radiation chemistry of water for TEM conditions by implementing a one-dimensional (1D) time-dependent reactiondiffusion model. The model comprises two types of species, namely primaries and secondaries, and describes the reaction kinetics in water upon continuous e --beam irradiation [13]. Primaries are defined as those species predominantly generated directly by irradiation (i.e., having non-zero generation cross-sections due to electron impact-G-values) [10], while secondaries are created through chemical reactions between primaries as well as other secondaries [13]. We refer the reader to other manuscripts for an in-depth description of both assumptions in the established radiolysis model (e.g., G-value theory and early heterogeneous stage), and its physical background [10,13,[22][23][24].
The initial model has been under continuous improvement and extension and has been adapted by several research groups aiming to better reflect the complex experimental realities. Existing implementations address aspects including the experimental parameters (e.g., temperature) [7,13,25], geometrical considerations [26], extending the chemical reaction system (beyond pure water) [22,[27][28][29][30], spatial inhomogeneities of energy deposition due to inhomogeneities of the scattering matter [22,26].
The development of general concepts for the mitigation of radiolytic effects is less advanced. One widely discussed concept involves additives, called scavengers, that are intentionally added to the reaction media and undergo fast reactions with the primary species decreasing their concentration [22,31,32]. By this means, the investigated process is expected to be less affected by radiochemistry-at least for a limited time span until the initial scavenger quantity is consumed. The main problem of this approach is that it is difficult to generalize since the peculiarities of each chemical system must be considered separately and appropriate scavengers must be found for every process individually [31,32].

Flow in Liquid-Phase TEM
Recently, a more general approach to mitigate radiolytic species in the irradiated volume was proposed in the LP-TEM community, namely: application of flow in the liquid flow cell (LFC) for continuous renewal of the reaction media during the imaging process [7,17,33,34]. To our knowledge, neither systematic experimental nor theoretical works were conducted so far to validate the approach. Thus, it remains unclear whether convective mass transport could serve for effective removal of the radiolytic species from the irradiated volume.
Various experimental LP-TEM setups supporting flow of the liquid in the imaging area are available nowadays [17,18,33,35,36]. All setups are based on microchips being allocated in the tip of a microfluidic TEM sample holder and can be categorised in two classes based on the way how the flow is directed. In the so-called 'bathtub' configuration, the microchips are immersed into an oversized flow channel and the liquid exchange in the LFC, and thus in the imaging area, is achieved mainly by diffusion [37,38]. In the second setup which is often referred to as 'nanofluidic', all the flow is forced into the nanochannel and liquid exchange in the imaging area is dominated by convection [35,38]. Recently it was shown that solutionexchange (due to diffusion) in a 'bathtub' setup is very slow and thus can hardly be utilised for liquid replenishment in the imaging area [38]; in an ideal 'nanofluidic' configuration, however, the linear flow velocity in the nanochannel can reach 1 μm/μs (according to information provided by the manufacturer) which is promising for efficient flushing of radiolytic species [18].
In this manuscript, we theoretically evaluate the influence of liquid flow in the irradiated area on the kinetics of the radiolytic reaction network and as a result on the concentration of reactive species. For this purpose, an existing radiolysis model [13] was implemented in a finite element (FE) simulation environment and extended by the convection components. In the range of experimentally reasonable velocities, the flow was found to have counterintuitive effects on the radiolytic reaction network in the irradiated volume with intriguing implications for LP-TEM experiments.

The Model
In the presented model, the well-established implementation for reaction-diffusion kinetics of radiolysis by Schneider et al. was extended with the physics of laminar flow in 2D space [13]. Equation 1 denotes the differential mathematical expression for calculating concentration fields, c i (x,y,t), of radiolytic species i in space, (x,y), and The right-hand terms in equation 1 correspond to (1) the volumetric creation rate due to radiolysis, (2) production and (3) consumption due to chemical reactions, (4) diffusivity of a species i and (5) convection due to continuous flow. The variables , r , y G i , F, k jk & k ij , D i and v denote the density of the irradiated liquid, absorbed dose rate, generation rate of species i due to absorbed dose (note, that G i = 0 for all but the 8 primary species; compare table S1 in Supporting Information), Faraday constant, rate constants of chemical reactions, diffusion coefficient of species i and the linear velocity of flowing liquid, respectively.  and 2  are 2D spatial gradient and Laplacian operators, respectively.

Geometrical description
The geometry of the model represents the nanochannel around the imaging area in the LP-TEM setup as a 2D cross-section in the plane along the flow direction and perpendicular to the beam path (refer to the Supporting Information section 3 and 6 for the justification of sufficiency of 2D simulations and related assumptions). All geometrical aspects of the model are defined in figure

Radiolysis and diffusion
The radiolytic creation rate of primary species (e h -,  H H . , OH  , H 2 , H 2 O 2 , H + , HO 2 . , OH − ; 1 st right-hand term in equation 1) within the irradiated region is calculated from constant non-zero G-values (A number of molecules created per 100 eV of absorbed energy, table S1) and an electron flux.
Upon creation, the species populate a complex chemical reaction network, that in neat water alone consists of 73 reactions and comprises 15 chemical species (8 primaries and 7 secondaries) and H 2 O [13]. The reactions are implemented as a set of kinetic equations (compare table S2) considering water as a solvent with much higher constant concentration [13].
Diffusive flux is an inherent consequence of concentration gradients due to the heterogeneous irradiation. Diffusion is mathematically described by Fick's 2 nd law (4 th term in equation 1).
The model was first implemented in 0 (no diffusion) and 1D (heterogeneous irradiation and diffusion) and tested against the reference model [13] (see Supporting Information section 4), and only then combined with the physics of the flow in a 2D channel geometry. Laminar flow Flow of a liquid through a channel (neglecting gravitational and inertial forces) results from a pressure gradient along the flow channel and depends strongly on the internal channel geometry and applied flow conditions [39][40][41][42]. The flow type, either turbulent or laminar, can be estimated by the characteristic Reynold number, which when being below ∼2100 will indicate a laminar flow [43]. In typical microfluidic scenarios including LP-TEM cells, the Reynolds number is well below 1 [38], indicating laminar flow (see Supporting Information section 1 for details).
In laminar flow, the liquid travels on parallel flow lines following a pressure gradient between in-and outlet; transport between flow lines occurs by diffusion only [41,44]. The flow velocity in the channel varies from the walls to the centre with maximum velocity at the centre and zero flow at the walls resulting from viscous forces (no-slip conditions) [39][40][41]. For laminar flow, the velocity profile between two infinitely long parallel plates separated by distance h is described by plane Poiseuille flow as [42] v y Q h where h and y are as indicated in figure 1. Q is a volumetric flow rate; v¯= Q/A is an average flow velocity derived from the cross-section A of the flow channel. Equation 2 describes a parabolic velocity profile with the maximum flow velocity, v max , in the middle of the channel equal to 3/2 of the average velocity, v. The relation between vā nd the experimental parameter Q can be established by 3D finite elements simulations based on the real geometry for every setup [38]. Since volumetric flow is obsolete for 2D geometry, we will focus on the more general parameter v¯at the inlet later in the text. Refer to Supporting Information section 1 for further details on fluid flow (in LP-TEM setups).

Boundary and initial conditions
The initial concentration of all species, except for H 2 O, H + and OH -, were zero in all domains, H + and OHconcentrations corresponded to pH and pOH 7. The production rate of primary species in the IA ( figure 1) was defined by corresponding G-values and an absorbed dose rate of 7.5 10 7 Gy s y = ⋅ (equivalent to 5.7 pA beam current assuming circular beam shape with the diameter W IA ).
The channel walls (black lines in figure 1) were considered impermeable and chemically inert (i.e., neglecting the effect of reactive surface sites, N s,i , on bulk solute concentration, c i = N i /V, which is reasonable since N s,i = N i ) for all species and no interaction with the electrons of the beam were considered (i.e., secondary electrons generated in the membranes were neglected which was demonstrated reasonable for materials with low electron scattering power such as SiN x ) [26]. No-slip conditions were assumed for the flow at the walls (i.e., the flow velocity was zero at the walls). The flow velocity was implemented as uniform mean inflow velocity v¯at the inlet (left edge in figure 1) directed into the channel. The outlet (right edge on figure 1) was set as a drain for both convective and diffusive mass transport. From available literature, a reasonable v¯inside the nanochannel of LP-TEM flow devices was estimated to be in the range from 0 to 1 m/s [18]. Refer to the Supporting Information section 1 and 3 for the mathematical implementation of the model and justification of the flow velocity range.
Solving the model The model was implemented in COMSOL Multiphysics ® software and solved with a finite element method [45]. A rectangular mesh was generated, and the optimal element size (2 nm × 2 nm spatial resolution) was determined in a convergence study (refer to Supporting Information section 4 for details).
In this manuscript, the stationary solution of equation 1 (i.e., dc i /dt = 0), which reflects the chemical equilibrium of all radiolytic species under the influence of constant irradiation in IA and fluid flow, is analysed. The stationary state is achieved in the time scales of ∼1 s according to the reference 1D model [13] and the equivalent time dependent solution for the 2D in-plane model obtained in this work (Supporting Information section 6 for details).

Results and discussion
The steady-state spatial distribution of all radiolytic species was calculated under constant irradiation in the IA and constant flow at mean velocities ranging from 0 to 1 m/s.  OHand OH .  H increase by 15%-30%. Changes in the concentration start at 10 −5 m/s and develop differently for different species with increasing velocity. All species experiencing significant increase of concentration are primaries, i.e., those with non-zero Gvalues. It is interesting to note that concentrations of H + and OHdo not change synchronously, pointing to the possibility of controlling acidity in the IA by flow velocity. As it was pointed out recently [46], H + and OHare out of equilibrium if water is irradiated, so radiolytic acidity π * should be used instead of pH. Dependence of this value on the flow velocity is provided in Supporting Information section 5.

Analysis of concentration increase
Assuming that the production rate of H . and e h in chemical reactions is relatively low as compared to radiolytic generation rate, the increase of concentration for the above discussed primary species is caused mostly by the decrease of the consumption rate. This decrease in consumption is the consequence of a reduced concentration of secondary reactants involved in the consumption reactions.
A potential reasoning to explain the increase of steady state concentration of certain species with an increase of flow velocity may be referred to as 'spatial reaction network separation'. In the presence of a flow, the spatiotemporal evolution of reaction networks is altered in respect to zero flow. Chemical reactions, which in the absence of the flow would occur at one and the same location inside the irradiated region (i.e., in the IA), are D the species has time to be rinsed from the IA and participates in the reaction network downstream; if τ i < t, D this species will react inside the IA. Lifetimes of the species may be defined as the time necessary to decrease its concentration to zero assuming zero order consumption reaction: cons is a sum of consumption rates for species i (= 3 rd right-hand term of equation 1). It should be noted that equation 3 is equivalent to the definition of the half-life time for the zero, first and second order reactions up to the constant coefficient ½ for zero and ln2 for the 1 st order reactions (see Supporting Information section 2 for detailed discussion) [47,48].  [13] follows that this steady state is established within units of seconds even for a 200 μm large model (see Supporting Information figure S21) and much faster for smaller closed volumes with limited outward-diffusion (steady state in the absence of diffusion occurs within ∼1 ms) [13]. It is interesting to compare qualitatively the diffusion fluxes of the species in the steady state at zero flow (see Supporting Information figure S10). Most of the species diffuse out of the IA, meaning that they are generated there in excess with respect to the surrounding. Only H 2 O 2 and O 2 diffuse back towards the IA, indicating that they are preferentially generated in excess outside the IA and are consumed inside the IA. The fastest consumption reactions for H 2  When applying flow, the concentration of accumulated species dramatically drops at flow velocities as small as 10 −5 m/s (see figure 3(a)), because even the minimal flow outperforms the diffusional mass transport, and the It is peculiar to note, that such a trend has been evidenced already by Schneider et al [13] (see e.g., figure 3 in referenced manuscript), who reported a high concentration of primaries at the timescales of 10 −5 s after irradiation started, which decreased later, when the concentration of secondaries (and thus consumption reaction rates) increased.
To validate this scenario, simulations for the flow of aerated water with dissolved O 2 at a concentration corresponding to the steady-state value of the initial simulation at zero flow (0.18 mM) were performed. Figure 5 shows the steady-state concentrations of O 2 and H . in the presence of such flow (compare to figure 2(b), (c); the colour range is the same on both figures). A straightforward consequence of the presence of O 2 in water is the almost complete suppression of the H . concentration increase. Moreover, an expected drop of O 2 concentration is observed downstream upon passing through the IA, confirming O 2 consumption in the IA. At 0.001 m/s diffusion still plays a role in the mass transport and is affecting the solution upstream.
It should be noted here that scavenging is a selective process. Thus, decrease in H . and e h concentrations by flowing O 2 -saturated water is accompanied by the increase of concentrations of H + , HO 2 and O 2 -(compare figure S17). In turn, application of H 2 O 2 flow (figure S14) suppresses more efficiently e h than H .  H as compared to O 2 (see figure S17 and S18).
General remarks Our numeric study shows that radiolysis is represented by a complex and mutable reaction network and its steady state strongly depends not only on the dose rate as was pointed out previously [13], but also on the mass transport in the system. In the analysis of distinct irradiation scenarios, not only the dose rate and reaction set should be considered, but also a complete geometrical representation of the experiment defining the spatial distribution of the reactants.
The simulations show that accumulation of such species as O 2 and H 2 O 2 in the liquid cell due to limited mass transport (e.g., static cell) leads to a self-scavenging of highly reactive atomic hydrogen and hydrated electrons. In this respect, the main action of the flow is the suppression of the accumulation of self-scavengers,