M G Pala
email: demetrio.logoteta@im2np.fr

Impact of the gate and insulator geometrical model on the static performance and variability of ultra-scaled silicon nanowire FETs Demetrio Logoteta, Nicolas Cavassilas, Alessandro Cresti, Marco G. Pala, Member IEEE, Marc Bescond

Abstract-We investigate the effect of the geometrical model adopted for the gate electrode and for the insulator envelopping the access regions on the full-quantum simulation of ultra-scaled nanowire field-effect transistors (NW-FETs). We compare the results obtained in the 'minimal' geometry commonly used in simulations with those obtained in a more realistic one, able to fully account for the gate fringing effects. We evaluate the impact of the model geometry on the static performance of NW-FETs, and discuss the interplay with the surface roughness and the random distribution of dopants. We find that the ION/IOFF ratio evaluated in the minimal geometry can be remarkably underestimated in short devices, notably in the case of small length-towidth ratio. The roughness-induced current degradation and the sensitivity to the surface roughness variability can also suffer from non-negligible underestimations when evaluated in this geometry. Finally, we point out that an inaccurate description of the device electrostatics is expected to result in an overestimation of the sensitivity to the doping-induced variability.

Index Terms-nanowires, fringe field, non-equilibrium Green's function, ION/IOFF ratio, surface roughness, doping induced variability.

I. INTRODUCTION

I NTEGRATED electronic devices and prototypes in exper- imental setups are commonly deeply buried into insulating layers and contacted by means of gate electrodes several tens of nanometers thick. Notwithstanding, the geometries routinely considered for the numerical quantum simulation of nanoelectronic devices assume for the insulator enrobbing the device access regions a quite small thickness with respect to the decay length of the electric field (usually the same thickness as the gate oxide) and model the gate with thickless stripes. This prevents any account for the outer fringe field originating from the gate sidewalls [START_REF] Kim | Modeling and significance of fringe capacitance in nonclassical cmos devices with gate-source/drain underlap[END_REF]. More realistic geometries are systematically considered only in the framework of simpler transport models, in order to provide semianalytical descriptions of parasitics [START_REF] Kim | Modeling and significance of fringe capacitance in nonclassical cmos devices with gate-source/drain underlap[END_REF], [START_REF] Liang | Fringe-induced barrier lowering (fibl) included sub-threshold swing model for double-gate mosfets[END_REF], [START_REF] Venugopalan | Bsim-cg: A compact model of cylindrical/surround gate mosfet for circuit simulations[END_REF] or perform extensive explorations of the design parameter space [START_REF] Sachid | Gate fringeinduced barrier lowering in underlap finfet structures and its optimization[END_REF], [START_REF] Zhao | Analysis of the effects of fringing electric field on finfet device performance and structural optimization using 3-d simulation[END_REF], [START_REF] Zhuge | Investigation of parasitic effects and design optimization in silicon nanowire mosfets for rf applications[END_REF], [START_REF] Wang | Investigation on variability in metal-gate si nanowire mosfets: Analysis of variation sources and experimental characterization[END_REF].

Here, we propose a full-quantum simulation study focusing on the static characteristics of gate-all-around ultra-scaled silicon NW-FETs, aimed to investigate the impact of extending the simulation domain so as to include a larger part of the gate electrode and of the insulator enveloping the source and drain extensions. Throughout the paper, we assess this impact in connection with some of the main tasks the numerical simulation of electronic devices is supposed to accomplish: the prediction of the upper performance bounds of devices, the estimation of their change in the presence of non-idealities, and the evaluation of their sensitivity to the sources of variability. We proceed by carrying out a systematic comparison of the results obtained for NW-FETs of different width and/or gate length within two different model geometries, which feature a different thickness of the gate and the insulator enrobbing the device. The deviations are quantified in terms of the I ON /I OFF ratio, which we assume as a representative figure of merit. Besides considering perfect devices, we also discuss the modifications of the I ON /I OFF ratio in the presence of surface roughness and its sensitivity to the disorder induced by surface roughness or random distribution of dopants. These latter are expected to be the dominant sources of variability for NW-FETs and have a stronger impact in ultra-scaled devices because of the weaker self-averaging effects.

The paper is organized as follows. Section II illustrates the physical models and the simulation techniques; Section III presents the main results: subsection III-A focuses on the evaluation of I ON in perfect NW-FETs; subsection III-B treats the current degradation induced by the surface roughness; in subsection III-C, we discuss the device sensitivity to the variability induced by surface roughness and doping. Finally, in Section IV we draw our concluding remarks.

II. MODEL

We consider square gate-all-around unstrained silicon NW-FETs along the [100] orientation. The dielectric constant of the gate insulator is set to 3.9 (SiO 2 ) and its thickness to t ox = 1 nm. This choice is aimed at simplifying as much as possible the device structure and make the interpretation of the simulation results easier. Replacing the SiO 2 layer with a more realistic, equivalent high-k stack is expected to negligibly change the device electrostatic, particularly the gate fringe field. However, a single SiO 2 layer of the considered thickness could entail a significant gate leakage current. This current could mask some of the effects emerging in the subthreshold regime of the devices and make difficult the extension of the resulting physical picture to NW-FETs with more opaque gate insulators. In light of this, the leakage current is not taken into account in our simulations.

The geometries compared in this study are sketched in Fig. 1. The minimal geometry A is the simplest and most commonly used one. The gate corresponds in it to a wrapping surface belonging to the boundary of the simulation domain. For a given size of the nanowire and gate oxide thickness, it is the rectangular parallelepiped geometry enclosing the minimum volume. On the other hand, the geometry B models the gate as a three-dimensional shell. The value of 5 nm set for its thickness has been verified to guarantee a reasonable stabilization of the results against any further increase. In both the A and B geometries, the spacer insulator is SiO 2 . This choice has been recently suggested as a better option to limit the increase of the switching delay as the gate length is scaled down to 10 nm [START_REF] Sachid | Bulk finfet with low-k spacers for continued scaling[END_REF]. In order to further simplify the device structure, SiO 2 was also set as the insulator enveloping the access regions of the devices, since the value of its dielectric constant is not far from the one of the commonly employed interlayer dielectrics.

In order to investigate the interplay with different degrees of quantum confinement and short channel effects, we considered devices with two different transversal cross-sections W 2 (2 × 2 and 5 × 5 nm 2 ), gate lengths L G (5 and 10 nm) and gate underlaps L S (0 and 3 nm). The source and drain extensions are assumed 12 nm long and continously n-doped with concentration N D = 10 20 cm -3 , unless otherwise stated. Following the projections of the International Technology Roadmap for Semiconductors 2015 for low power devices at the 2021/2024 technological nodes [START_REF]International technology roadmap for semiconductors[END_REF], throughout the paper we assume a source-drain bias of V DS = 0.6 V and we set the nominal OFF current to I OFF = 10 -4 µA/µm, unless otherwise stated. The current in the ON state is computed as

I ON = I DS (V (off) GS + V DS ), where V (off) GS
denotes the gate voltage V GS in the OFF state. All the values of current provided in the paper are normalized to the width of the nanowires.

Our simulation methodology is based on the self-consistent solution through a fixed point iteration algorithm of the transport equations within the non-equilibrium Green's function formalism and of the Poisson equation. A comprehensive illustration of this approach can be found in Ref. [START_REF] Pourfath | The Non-Equilibrium Green's Function Method for Nanoscale Device Simulation[END_REF]. To describe the system, we assumed an effective mass Hamiltonian, calibrated on tight-binding band structure computations [START_REF] Nehari | Influence of band structure on electron ballistic transport in silicon nanowire mosfets: An atomistic study[END_REF]. This approximation is able to provide an accurate description of transport in silicon nanowires with the considered crystallographic orientation [START_REF] Neophytou | Simulations of nanowire transistors: Atomistic vs. effective mass models[END_REF]. The nanowires with a 2 × 2 nm 2 transversal cross-section have been simulated in the realspace basis corresponding to a seven-point stencil central finite difference discretization. For the nanowires with W = 5 nm a coupled mode-space approach has been resorted, in order to limit the otherwise unmanageable computational burden. Electron-phonon scattering was included within the selfconsistent Born approximation. Acoustic and optical phonons were modeled in the elastic and dispersionless approximation, respectively, by assuming the self-energy expressions and the deformation potentials in Ref. [START_REF] Bescond | Modeling of phonon scattering in n-type nanowire transistors using one-shot analytic continuation technique[END_REF].

The surface roughness at the interface between silicon and SiO 2 was included nonperturbatively by introducing random fluctuations in the oxide thickness, according to the metodology detailed in Ref. [START_REF] Buran | Three-dimensional real-space simulation of surface roughness in silicon nanowire fets[END_REF]. The surface fluctuations ∆(r) have been generated by assuming an exponential autocorrelation function

∆(r)∆(r ′ -r) = ∆ 2 m e - √ 2r/Lm
, where ∆ m is the root mean square (rms) value of the fluctuations and L m their correlation length.

Ionized discrete donors have been modeled as positive elementary point charges. Central cell corrections, depending on the particular chemical species considered and expected to be small for shallow dopants [START_REF] Yu | Fundamentals of Semiconductors[END_REF], have been neglected. In order to take into account the substitutional position of dopants, a hard sphere model has been resorted, where individual dopants have been considered at the center of nonoverlapping spheres with radius 1 nm. The screening by the mobile charge is automatically taken into account within the Hartree approximation through the self-consistent coupling with the Poisson equation.

L G =5 nm, A L G =5 nm, B L G =10 nm, A L G =10 nm, B L G =10 nm BAL, A L G =10 nm BAL, B 0 0.2 0.4 V GS (V)
V DS =0.6 V W=5 nm V DS =0.6 V L S =3 nm L S =3 nm x 10

III. RESULTS AND DISCUSSION

A. ON current in perfect NW-FETs

Fig. 2 and Fig. 3 show the transfer characteristics for the NW-FETs with L S = 3 nm and L S = 0, respectively, in the A and B geometry. All the curves are shifted in order to set I DS = I OFF at V GS = 0 V. Table I summarizes the percentage increase of the I ON /I OFF ratio in B with respect to A. It turns out that the I ON is always higher in B, by an amount depending on the particular configuration considered. The largest variations are obtained for the devices with W = 5 nm and L G = 5 nm, for which the ON current can Energy (eV) double (L S = 3 nm) or even increase several times (L S = 0). According to Fig. 2 and Fig. 3, this behavior is due to the significant subthreshold swing degradation in A with respect to B. More insight can be gained from Fig. 4, showing the lowest conduction subband in the NW-FET with W =5 nm, L G = 5 nm and L S = 3 nm. The plot highlights that for a given V GS the source-to-drain barrier is wider and higher in B, which can be interpreted as due to an increase of the effective gate length with respect to A. This behavior can be traced back to the different extent by which the geometries account for the gate fringe field. Geometry B is able to fully account for it (inner fringe and outer fringe from the gate sidewalls [START_REF] Kim | Modeling and significance of fringe capacitance in nonclassical cmos devices with gate-source/drain underlap[END_REF]), while geometry A can only take into account the inner fringe. This is further illustrated in Fig. 4 by the plot of the difference between the subbands obtained in A and B corresponding to the same V GS : the difference peaks in the source and drain spacers, where the fringe field is more intense, and decreases as V GS increases, due to the enhancement of the gate screening with the increment of the mobile charge. Analogously to the case in which the effective gate length is intentionally increased through the inclusion of gate underlaps [START_REF] Trivedi | Nanoscale finfets with gate-source/drain underlap[END_REF], this effect entails an increase of the control over the height of the barrier and a suppression of the tunneling current in the subthreshold region, thus resulting in an improved subthreshold swing in the B geometry. The relative improvement is higher for the devices with larger W , shorter gate length and no spacers, since they exhibit a poorer electrostatic integrity and the relative impact of an increase of the effective gate length is stronger. We point out that the discrepancies between the value of the I ON /I OFF ratio evaluated in the A and B geometry can strongly depend on the electron-phonon coupling strength. The percentage increase of the I ON /I OFF ratio in the ballistic regime is reported in Table I for the NW-FETs with L S = 3 nm. The comparison with the values obtained in the presence of electron-phonon scattering indicates that the deviations between A and B are always larger in the ballistic regime, in particular for the narrower NW-FETs. As can be checked in Table I by comparing the cases L G = 5 nm and 10 nm, the magnitude of the effect is weakly dependent on the gate length, which suggests that it can be traced back to a different modulation of the potential at large gate overdrives. In order to gain more insights, the lowest conduction subbands on the source side corresponding to the transfer characteristic data points enclosed by a rectangle in Fig. 2, are compared in Fig. 5. In both the ballistic and dissipative case, switching from A to B geometry entails a decrease of the potential barrier in the source spacer, which is directly modulated by the fringe field. On the other hand, while in the ballistic case electrons are mainly injected from the source contact at the base of the barrier, in the presence of phonons electrons are injected over its top, due to the linear bending of the subband (see maps in Fig. 5). Only a fraction of them hit the barrier due to the incomplete energy relaxation in the source. As a results, the modification of the barrier in B with respect to A at large gate overdrives affects more the ballistic current than its dissipative counterpart. A similar behavior is found for the NW-FETs with W = 2 nm and L S = 0 (not shown). In wider devices, the modulation of the barrier at large gate overdrives is less influenced by the gate fringe field, and similar percentage increases of the I ON /I OFF ratio in B with respect to A are found in the ballistic regime and in the presence of electron-phonon scattering (cf. the values for W = 5 nm in Table I). This is a consequence of the lower resistance of the source extension, which results in a smaller bending of the subbands, and of the weaker effectiveness of the fringe field in modulating the potential over the whole transversal cross-section of the nanowire.

L G =5 nm, A L G =5 nm, B L G =10 nm, A L G =10 nm, B W=2 nm L S =0 V DS =0.6 V x 10 2 (a)
B A B-A W=5 nm L G =5 nm V GS Gate Spacer Spacer L S =3 nm

B. Surface roughness-induced current degradation

We limit our analysis to NW-FETs with W = 5 nm. In NW-FETs with W = 2 nm typical rms values for the surface roughness result in considerably larger subband fluctuations, entailing a stronger dependence of the transfer characteristics on the particular roughness realization. In these conditions a full statistical study is required in order to extract clear trends. Fig. 6 illustrates the percentage degradation of the current as a function of V GS for NW-FETs with L G = 5 nm and L S = 3 nm and for rms values of 0.2 and 0.4 nm. In order to make the comparison meaningful, the same realization of the surface roughness has been used for A and B. The plots clearly show a larger current degradation in B, which reduces as V GS increases. We discuss first the low-V GS behavior. At low gate voltages only few low-energy modes contribute to transport [START_REF] Poli | Size dependence of surface-roughness-limited mobility in silicon-nanowire fets[END_REF]. Due to the presence of a large potential barrier in the channel, in the source extension the transport via these modes takes place at higher energies with respect to the bottom of the corresponding subbands and is thus only weakly affected by the subband fluctuations induced by the surface roughness. The roughness influences the transport mainly in the channel region, by affecting the shape and the height of the barrier. On the other hand, as already mentioned in Sec. III-A and quantified in the inset of Fig. 6, the source-drain tunneling is overestimated in the A geometry. The tunnel component of the current at low enough energies is weakly sensitive to the details of the top of the barrier and thus, in the considered situation, to the fluctuations induced by the roughness. Therefore, the overestimation of tunneling finally results in an underestimation of the scattering and, as a consequence, of the current degradation in A with respect to B. The weak sensitivity of the low-energy tunneling current to the potential fluctuations can be understood by referring to the Wentzel-Kramers-Brillouin (WKB) approximation for the transmission through a barrier defined by a potential U (x) = U 0 (x) + δU (x), where U 0 is the smooth potential and δU (x) the potential fluctuations due to the roughness. The transmission T as a function of the energy E reads T (E) ∼ e -(2/ ) √ 2m * (E-U (x))dx , where is the reduced Plank constant, m * is the effective mass in the transport direction and the integral is computed between the classical turning points. To the first order in the fluctuations δU (x), the variation of the action integral is δ

2m * (U (x) -E)dx ∝ δU (x) √ (U0(x)-E)
dx. This shows 0 0.2 0.4 0.6

V GS (V) that the influence of a perturbation δU (x) on the transmission tends to become negligible as E decreases with respect to the top of the barrier.

As V GS increases, the discrepancies between the A and B case reduce. This originates from the concomintant decrease of the percentage of tunneling and of the difference of this percentage between A and B (see the inset of Fig. 6). Moreover, as the components of the current spectrum close to the bottom of the subbands increase due to the lowering of the barrier, also the surface roughness in the source and drain extension starts to contribute to the back scattering, thus reducing the relative importance of the fluctuations in the channel. The deviations in the current degradation entails the emergence of deviations in terms of the I ON /I OFF ratio. Since, according to the previous discussion, I OFF in rough devices is degraded more than I ON , it turns out that the I ON /I OFF ratio is enhanced by the surface roughness (both in A and B). For ∆ m = 0.4 nm, we find for the relative variation of the I ON /I OFF ratio with respect to the smooth case |∆(I ON /I OFF )|/(I ON /I OFF ) smooth = 72% in A, while it appreciably increases to 96% in B.

For longer gates (not shown), the difference in the current degradation between A and B decreases due to the suppression of the tunneling.

C. Device sensitivity to the sources of variability

In this section we consider the modifications of the NW-FET transfer characteristics induced by a change of the surface roughness realization or of the dopant distribution. The purpose is not to assess the variability ranges, which would require a full statistical study, but only to qualitatively investigate how the sensitivity predicted for the devices changes when switching from the A to the B geometry.

1) Surface roughness-induced variability: In order to address the impact of the model geometry on the sensitivity to the surface roughness variability, an ensemble {1, 2} of two different realizations (the same for A and B) has been considered. The relative deviation of the current ξ = |I

DS + I (1)
DS )/2 is taken as a measure of the sensitivity. ξ is plotted in Fig. 7 for NW-FETs with W = 5 nm, L G = 5 nm and L S = 3 nm. A V GS dependence qualitatively similar to the one associated to the current degradation in Fig. 6 is found, with a clear enhancement of ξ in B with respect to A. This translates to a different sensitivity of the OFF current and thus of the I ON /I OFF ratio. By defining the relative deviation of the I ON /I OFF ratio with respect to the average as

ξ ION/IOFF = I (2) ON I (2) OFF - I (1) ON

I

(1) OFF

I ON I OFF , (1)
we find ξ

(A)
ION/IOFF = 2.6% and ξ (A)

ION/IOFF = 11%, for rms=0.2 nm and 0.4 nm, respectively. The corresponding values in B, 6.4% and 20%, increase by almost a factor 2.

As for the current degradation in the previous subsection, the origin of the discrepancies in Fig. 7 can be traced back to the overestimation of the tunneling in the geometry A. Indeed, in the previous subsection we have pointed out that the roughness has a weaker impact on the current in A. This directly entails a smaller sensitivity of the current the change of the roughness realization with respect to B. On the other hand, by referring again to the WKB approximation for the tunneling transmission, we can write

T (E) ∼ e -(2/ ) √ 2m * (U (x)-E)dx = e -( √ 8m * / )∆x √ (U (x)-E)
x , where ∆x is the width of the barrier at the energy E and f (x) x denotes the average of f (x) over ∆x. The rearranged expression for T indicates that the tunneling current tends to experience the potential in an averaged way, or, in the examined situation, that the tunneling current can be approximately considered modulated by the potential averaged over the roughness fluctuations, instead than by the actual rough potential. This further contributes to decrease the current sensitivity to the roughness in A with respect to B. The dependence of the sensitivity on V GS in Fig. 7 is qualitatively similar to the one of the current degradation in Fig. 6, and the same considerations made in the previous subsection apply.

2) Doping-induced variability: We consider NW-FETs with L S = 0 and the structure sketched in Fig. 8, including on the source and drain side a 6 nm long region where the continous doping is replaced with a distribution of positive point charges modeling ionized dopants. The number of charges has been set equal to the average number of discrete dopants in the considered volume region according to the doping concentration N D . While this approach entails an unphysical interface between the continuously doped region and the one including discrete dopants, it guarantees that the (screened) Coulombic tails of the dopant potential extending toward the source and drain contacts are fully taken into account. An alternative approach to model disordered contacts can be found in Ref. [START_REF] He | Non-equilibrium green's functions method: Non-trivial and disordered leads[END_REF]. We considered two representative arrangements of dopants, differing in the average distance of the dopants from the gated region. Since the dopants close to the channel are very effective in lowering and thinning the source-todrain barrier, this distance is expected to be one of the dominant parameters determining the dispersion of the transfer characteristics [START_REF] Martinez | Quantum-transport study on the impact of channel length and cross sections on variability induced by random discrete dopants in narrow gate-all-around silicon nanowire transistors[END_REF]. The transverse y and z coordinate of each The empty and filled circles of the plot correspond to the projection on the x -y plane of the two distributions of discrete dopants considered in the source extension of the NW-FET with W = 5 nm. Empty circles have a Gaussian probability density of their x coordinate (dashed curve, not normalized) centered at a distance δ = 4 nm from the gated region; in the case of the filled circles the probability density (continue curve) is shifted of 2 nm in the direction of the gate. The y and z coordinates of the dopants are extracted for both distributions according to a uniform probability density. point charge have been independently extracted assuming a uniform probability distribution over the silicon transversal cross-section. The x coordinate of the charges follows instead a Gaussian distribution. The two different arrangements of dopants correspond to different choices for the mean of the Gaussian probability density, namely for the x coordinate x bar of the barycenter of the charge distribution. For one of the arrangement, which we will denote by 'F' (Far), x bar is 4 nm far from the gated region, both on the source and the drain side, while it is 2 nm far for the other, which we will denote by 'C' (Close). Analogously to the case of rough devices in the previous subsection, we computed the transfer characteristics with the two dopant arrangements both in the A and the B geometry and we then extracted the relative deviation of the current with respect to the average ξ = |I Note that here we set a larger value of I OFF because for the considered dopant arrangements the deviations between A and B emerge, in terms of ξ, only at large enough gate overdrives. Indeed, in the devices with the dopant distributions closer to the gate, the subthreshold swing is strongly degraded and, in the subthreshold region, I Fig. 9 shows the plot of ξ in the geometry A and B for NW-FETs with W = 5 nm and W = 2 nm. The sensitivity appears always overestimated in A. This discrepancy originates from the antagonistic effect of the gate fringe fields, partially counterbalancing the differences in the electrostatic potential due to the different distribution of the dopants. Indeed, while the positive dopants close to the gate tend to thin and lower the potential barrier, the gate fringe field tends to widen and lift it. This interplay finally results in a weaker sensitivity of the self-consistent potential in the channel to the distribution of dopants, translating to a smaller sensitivity of the transfer characteristics. Since the gate fringe field gradually decreases far from the gate just as the effect of a dopant on the barrier decreases as it is moved apart toward the source or drain contact, this behavior can be interpreted as a disorder regularization promoted by the fringe field and tending to level out the influence of the dopants at different positions. Finally, we report the deviations in terms of ξ ION/IOFF , for which a definition analog to Eq. 1 is assumed. Despite the differences in the current at large gate overdrives visible in Fig. 9, no sensible deviations are found for the wider NW-FETs because, since I

IV. CONCLUSION

We compared the results of full-quantum simulations of ultra-scaled silicon NW-FETs obtained in a geometry fully accounting for the gate fringe field with the corresponding obtained in a commonly used minimal geometry only partially accounting for them. The comparison was carried out in terms of the prediction of the upper bounds of the device static performances, of their change in the presence of surface roughness, and of their sensitivity to the variability induced by surface roughness and random distribution of dopants. Our results indicate that large underestimations of the I ON /I OFF ratio can arise in the minimal geometry, especially for devices with high W/L G . The poor electrostatic integrity of these devices enhances, indeed, the relative importance of the effective gate length increase promoted by the outer fringe field. In conditions of weaker electron-phonon scattering (lower temperature or materials exhibiting a weaker electron-phonon coupling), this underestimation is expected to increase for very narrow NW-FETs due to the stronger sensitivity to the modulation of the barrier exerted by the fringe field at large gate overdrives. The current degradation and the sensitivity to the variability induced by the surface roughness at the channel-oxide interface are predicted to be non-negligibly underestimated in the minimal geometry. This effect has been traced back to the overestimation of the tunnel component of the current, less affected by the scattering from the fluctuations of the electrostatic potential. Finally, we highlighted a fringe field-driven regularization effect on the disorder induced by the random distribution of dopants, entailing an overestimation of the sensitivity to the doping-induced variability when the fringe field is not fully taken into account.

This study helps clarifying the physical mechanisms governing the interplay between the gate fringe fields and the electron transport, and provides a guide to estimate the degree of approximation associated to different choices of the simulated geometry.

Fig. 1 .

 1 Fig. 1. Sketch of the simulated NW-FETs, modeled according to the minimal geometry A and the more realistic geometry B.

Fig. 2 .

 2 Fig. 2. Transfer characteristics of the NW-FETs with L S = 3 nm, simulated with the A and the B geometry. (a) W = 2 nm; (b) W = 5 nm. All the curves have been shifted in order to set I DS = I OFF at V GS = 0 V. The abbreviation 'BAL' in the legends indicates that the corresponding curves refer to the ballistic regime.

Fig. 3 .

 3 Fig. 3. Transfer characteristics of the NW-FETs with L S = 0 nm, simulated with the A and the B geometry. (a) W = 2 nm; (b) W = 5 nm. All the curves have been shifted in order to set I DS = I OFF at V GS = 0 V.

Fig. 4 .

 4 Fig.[START_REF] Sachid | Gate fringeinduced barrier lowering in underlap finfet structures and its optimization[END_REF]. Lowest subbands of the NW-FET with W = 5 nm, L G = 10 nm and L S = 3 nm for three different gate voltages, simulated by using the A and the B geometry (dashed and dot-dashed lines, respectively). The continous lines represent the difference between two subbands in the A and the B geometry corresponding to the same V GS value.

Fig. 5 . 2 0Fig. 6 .

 526 Fig. 5. Current spectrum on the source side for NW-FETs with W = 2 nm, L G = 10 nm and L S = 3 nm simulated with the B geometry in (a) the ballistic regime and (b) in the presence of electron-phonon scattering.The continuous curves represent the corresponding lowest subband, while the dashed ones refer to the lowest subband obtained in the A geometry. The maps and the curves refer to the points enclosed by a rectangle in Fig.2

Fig. 7 .

 7 Fig. 7. Deviation of the current with respect to the average ξ = |I (2) DS -I (1) DS |/ I DS for a set {1, 2} of two different realizations of the surface roughness. Two different rms values of the oxide thickness fluctuations ∆m =0.2 and 0.4 nm are considered. Lm = 1 nm. W = 5 nm, L G = 5 nm and L S = 3 nm.

  |/ I DS with respect to the average I DS = (I (2)

Fig. 8 .

 8 Fig.8. Structure of the nanowires considered for the study of the doping variability. The empty and filled circles of the plot correspond to the projection on the x -y plane of the two distributions of discrete dopants considered in the source extension of the NW-FET with W = 5 nm. Empty circles have a Gaussian probability density of their x coordinate (dashed curve, not normalized) centered at a distance δ = 4 nm from the gated region; in the case of the filled circles the probability density (continue curve) is shifted of 2 nm in the direction of the gate. The y and z coordinates of the dopants are extracted for both distributions according to a uniform probability density.

  DS |/ I DS as a function of V GS . For transfer characteristics obtained in the same geometry a common reference is set by shifting them so that (I 2 = I OFF = 10 -1 µA/µm at V GS = 0.

  ξ (A) ≈ ξ (B) ≈ 200%.

Fig. 9 .

 9 Fig. 9. Deviation of the current with respect to the average ξ = |I (C) DS -I (F) DS |/ I DS for a set {C, F} of two different distributions of discrete dopants. W = 2 and 5 nm, L G = 5 nm.

  ≈ 200%. For the narrower NW-FETs the degradation of the subthreshold swing in the 'C' case is weaker and we find ξ (A) ION/IOFF ≈ 116% and ξ (B) ION/IOFF ≈ 126%.

TABLE I PERCENTAGE

 I INCREASE OF THE I ON /I OFF RATIO IN THE B GEOMETRY WITH RESPECT TO THE A ONE IN THE BALLISTIC (BAL) AND DISSIPATIVE (PHON) REGIME. I OFF IS SET TO 10 -4 µA/µm.

			LS = 3 nm		LS = 0 nm
		LG = 5 nm	LG = 10 nm	LG = 5 nm LG = 10 nm
		BAL	PHON BAL PHON	PHON	PHON
	W = 2 nm	53%	28%	26%	10%	26%	6%
	W = 5 nm 122%	117%	21%	20%	516%	13%

This work was funded by the French Agence Nationale de la Recherche via the project "NOODLES" (n. ANR-13-NANO-0009-01). D. Logoteta, N. Cavassilas and M. Bescond are with IM2NP


Dr Goedele Roos
email: goedele.roos@univ-lille.fr

Prof. Dr Jeremy N Harvey
email: jeremy.harvey@kuleuven.betoc

Histidine versus cysteine bearing heme-dependent halogen peroxidases: parallels and differences for Cl ➖ oxidation

The homodimeric myeloperoxidase (MPO) features a histidine as proximal ligand and a sulfonium linkage covalently attaching the heme porphyrin ring to the protein. MPO is able to catalyze Cl -oxidation with about the same efficiency as chloroperoxidase at pH 7.0. In this study, we seek to explore the parallels and differences between the histidine and cysteine hemedependent halogen peroxidases. Transition states, reaction barriers and relevant thermodynamic properties are calculated on protein models. Together with electronic structure calculations it gives an overview of the reaction mechanisms, and of the factors that determine the selectivity between one-and two-electron paths. Conclusions point to the innate oxidizing nature of MPO with the ester and sulfonium linkages hiking up the reactivity to enable chloride oxidation. The installation of a deprotonated imidazolate as proximal ligand shifts the does not equilibrium from one-to two-electron events without influencing the chemistry of the oxidation reaction.

Introduction

Heme peroxidases utilize hydrogen peroxide to oxidize organic compounds as an essential step in cell component synthesis and in metabolic and elimination pathways. They also form an essential part of the front-line defense against pathogens of the innate immune system by releasing highly reactive and oxidizing antimicrobial reaction products as hypohalides [START_REF] Klebanoff | Myeloperoxidase: Friend and Foe[END_REF] . In addition, MPO is increasingly reported to play a role in disease pathogenesis, because of its ability to generate reaction products that can oxidize many types of biomolecules such as lipids, DNA and proteins. For an overview of relevant literature, we refer to references given in Soubhye et al. [START_REF] Soubhye | Structure-Based Design, Synthesis, and Pharmacological Evaluation of 3-(Aminoalkyl)-5-Fluoroindoles as Myeloperoxidase Inhibitors[END_REF] and Aldib et al. [START_REF] Aldib | Evaluation of New Scaffolds of Myeloperoxidase Inhibitors by Rational Design Combined with High-Throughput Virtual Screening[END_REF] .

Peroxidases catalyze one-(Eq. ( 1)) and two-(Eq. ( 2)) electron oxidation reactions, corresponding to the so-called "one-electron" and "two-electron" pathways respectively [START_REF] Furtmuller | Active Site Structure and Catalytic Mechanisms of Human Peroxidases[END_REF] . H2O2 + 2AH2 ® H2O + 2AH i (1) H2O2 + X -+ H + ® H2O + HOX (2) In both pathways, the peroxidases undergo a series of redox reactions, hereby cycling in the two-electron pathway through the native enzyme (ferric MPO, Fe (III)) and compound I (oxoiron Fe(IV)=O porphyrin radical cation, Cpd I). In the one-electron pathway, Cpd I is recycled back to the native enzyme through compound II (protonated oxoiron Fe(IV)) (Figure 1). The central redox intermediate Cpd I, is the product of the two-electron oxidation of the native ferric heme by hydrogen peroxide (Eq. ( 3)).

Por ××× Fe(III) é ë ù û + H2O2® Our focus in this study is to explore the mechanism of two-electron chemistry involved in the oxidation of halides X ➖ (X= Cl, Br, I) to hypohalides (HOX) catalyzed by Cpd I (Eq. ( 4)), which is reduced during the reaction to the native enzyme (cfr. Eq. ( 2)). We also wish to explore which features make a peroxidase well-suited to carrying out this reaction, and determine the competition between the pathways.

i+ Por ××× Fe(IV)=O é ë ù û + X -+ H + ® Por ××× Fe(III)

é ë ù û + HOX (4)
Two electron chemistry has thus far mainly been reported for heme proteins possessing a negatively charged thiolate or tyrosinate as proximal ligand [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF] . The electron donating ligands elevate the pKa of Cpd II and thus lower the one-electron potential of Cpd I [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF][START_REF] Yosca | Iron(IV)Hydroxide PK(a) and the Role of Thiolate Ligation in C-H Bond Activation by Cytochrome P450[END_REF][START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] , hereby suppressing one-electron oxidation by shifting the equilibrium from sequential one-electron to two-electron events. This effect, also called the "pushing" effect of cysteine thiolate, by which futile one-electron cycles are avoided, was proposed as underlying mechanism by which C-H bond activation by cytochrome P450 is possible [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF][START_REF] Yosca | Iron(IV)Hydroxide PK(a) and the Role of Thiolate Ligation in C-H Bond Activation by Cytochrome P450[END_REF] and by which catalase can effectively disproportionate H2O2 into oxygen and water [START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] . For the same reasons, peroxidases with a (neutral) histidine ligand, having a low Cpd II pKa, typically operate through one-electron reactions [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF] . Nevertheless, a number of secreted chordata peroxidases such as myeloperoxidase (MPO), lactoperoxidase (LPO) and eosinophil peroxidase (EPO), which have histidine ligands, are also reported to carry out halogen oxidation. In this study, we seek to explore the parallels and differences between the histidine and cysteine heme-dependent halogen peroxidases.

CPO, a highly versatile heme peroxidase from the marine fungus Caldariomyces fumago, has a negatively charged cysteine-derived proximal ligand, similar to the related monooxygenase cytochrome P450 and to nitric oxide synthase (NOS). For an overview of literature on CPO, we refer to references given in Kühnel et al. [START_REF] Kühnel | Crystal Structures of Chloroperoxidase with Its Bound Substrates and Complexed with Formate, Acetate, and Nitrate[END_REF] . Halogen peroxidases oxidize X -to HOX with an apparent second order rate constant of 2.5 10 4 M -1 s -1 (Table 1), but are also able to catalyze many P450-like oxygen insertion reactions. In general, apparent second order rate constants for reactions of Cpd I with two-electron reductants catalyzed by CPO show attenuated activity compared to P450. For an overview, see Zhang et al. [START_REF] Zhang | Kinetics of Two-Electron Oxidations by the Compound I Derivative of Chloroperoxidase, a Model for Cytochrome P450 Oxidants[END_REF] .

The secreted chordata peroxidases 10 such as the homodimeric myeloperoxidase (MPO), the momomeric lactoperoxidase (LPO) and eosinophil peroxidase (EPO) feature a histidine as proximal ligand and are therefore unique in their ability to efficiently oxidize halides to hypohalides in a two-electron reaction (Table 1) [START_REF] Furtmuller | Active Site Structure and Catalytic Mechanisms of Human Peroxidases[END_REF] . The proximal histidine in these proteins has been suggested to be deprotonated (and thus negatively charged) [START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF] . MPO, EPO and LPO are all able to oxidize Br ➖ , I ➖ and SCN ➖ , but only the MPO clade is able to catalyze the production of hypochlorous acid (HOCl) from hydrogen peroxide (H2O2) and chloride anions (Cl ➖ ) with about the same efficiency (2.5 10 4 M -1 s -1 ) [START_REF] Furtmuller | Active Site Structure and Catalytic Mechanisms of Human Peroxidases[END_REF] as reported for CPO at pH 7.0 (Table 1). The HOCl/Cl -redox couple has an extremely high standard reduction potential E °' at pH 7.0 of 1.23 V, while those of other HOX/X ➖ redox couples are substantially lower: X=Br: 1.13 V; X= I: 0.78 V; X= SCN: 0.56 V. The rate at which Cpd I of MPO is able to oxidize these (pseudo-)halides follows the sequence of their reduction potential. 15 in Zhang et al. [START_REF] Zhang | Kinetics of Two-Electron Oxidations by the Compound I Derivative of Chloroperoxidase, a Model for Cytochrome P450 Oxidants[END_REF] CPO, MPO, LPO and EPO are heme b -containing peroxidase having a protoporhyrin IX as basic structure of its prostatic group. MPO's porphyrin ring is covalently attached to the enzyme via two ester bonds (Asp94 and Glu 242) and one electron withdrawing sulfonium linkage (Met 243) (Figure 2) [START_REF] Fiedler | X-Ray Crystal Structure and Characterization of Halide-Binding Sites of Human Myeloperoxidase at 1.8 A Resolution[END_REF] . LPO and EPO only feature the two ester linkages, while the porphyrin in CPO is not linked to the protein. This makes MPO featuring unique post-translational modifications of its heme group.

Figure 2 Protoporhyrin IX of MPO in its native ferric Fe(III) state and its covalent attachment to the protein via two ester bonds (Asp94 and Glu242) and one sulfonium linkage (Met243). In addition, the proximal His336 is shown.

Post-translational modifications impose unusual properties of the heme porphyrin. Resonance Raman spectroscopy indicates the enhancement of several out-of-plane low frequency modes by the presence of the ester bonds [START_REF] Zederbauer | Heme to Protein Linkages in Mammalian Peroxidases: Impact on Spectroscopic, Redox and Catalytic Properties[END_REF][START_REF] Brogioni | Resonance Raman Assignment of Myeloperoxidase and the Selected Mutants Asp94Val and Met243Thr. Effect of the Heme Distortion[END_REF] . The UV-Vis absorption spectra of MPO in the Fe(III) in the resting state are characterized by an intense Soret-peak at 428-430 nm [START_REF] Furtmuller | Active Site Structure and Catalytic Mechanisms of Human Peroxidases[END_REF][START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF][START_REF] Harrison | Compound I of Myeloperoxidase[END_REF][START_REF] Kooter | Site-Directed Mutagenesis of Met243, a Residue of Myeloperoxidase Involved in Binding of the Prosthetic Group[END_REF][START_REF] Kooter | The Sulfonium Ion Linkage in Myeloperoxidase: DIRECT SPECTROSCOPIC DETECTION BY ISOTOPIC LABELING AND EFFECT OF MUTATION[END_REF] , which is strongly red-shifted compared to the spectra of other mammalian heme b peroxidases having a proximal histidine ligand as LPO [START_REF] Furtmüller | Reaction of Lactoperoxidase Compound I with Halides and Thiocyanate[END_REF] and EPO [START_REF] Furtmüller | Spectral and Kinetic Studies on the Formation of Eosinophil Peroxidase Compound I and Its Reaction with Halides and Thiocyanate †[END_REF] (412-415 nm) and of the horseradish peroxidase ( ∼ 403 nm) 10,22 , having no post-translational modifications onto its protoporfyrine IX prostatic group (See Table 2 in Nicolussi et al. [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF] for an extensive overview of electronic absorption maxima of MPO, its mutants and several other peroxidases). This red-shift could be attributed by mutagenesis studies [START_REF] Zederbauer | Heme to Protein Linkages in Mammalian Peroxidases: Impact on Spectroscopic, Redox and Catalytic Properties[END_REF][START_REF] Brogioni | Resonance Raman Assignment of Myeloperoxidase and the Selected Mutants Asp94Val and Met243Thr. Effect of the Heme Distortion[END_REF][START_REF] Kooter | The Sulfonium Ion Linkage in Myeloperoxidase: DIRECT SPECTROSCOPIC DETECTION BY ISOTOPIC LABELING AND EFFECT OF MUTATION[END_REF] and by TD-DFT calculations of absorption spectra [START_REF] Devarajan | Effect of Covalent Links on the Structure, Spectra, and Redox Properties of Myeloperoxidase -A Density Functional Study[END_REF] to the presence of both the esters and sulfonium linkages. Further, post-translational modifications influence the redox properties of the heme. The presence of the positively charged sulfonium linkage has been linked to the high reduction potentials found for all relevant redox couples of MPO, compared to LPO and other heme b peroxidases (Table 2). For example, typically, the standard reduction potential E °' at pH 7.0 of the Fe(III)/Fe(II) redox couple of heme b peroxidases is negative (e. g. LPO [START_REF] Furtmuller | Standard Reduction Potentials of All Couples of the Peroxidase Cycle of Lactoperoxidase[END_REF] : -176 mV; HRP [START_REF] Battistuzzi | Redox Thermodynamics of the Fe 3+ /Fe 2+ Couple in Horseradish Peroxidase and Its Cyanide Complex[END_REF] : -306 mV) whereas this of MPO is positive (+5 mV) [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF][START_REF] Arnhold | Redox Properties of the Couple Compound I/Native Enzyme of Myeloperoxidase and Eosinophil Peroxidase[END_REF] . Reduction potential calculations of the FeII/FeIII and Cpd I/Cpd II redox couples show indeed that especially the sulfonium linkage shifts the potentials to higher positive values [START_REF] Devarajan | Effect of Covalent Links on the Structure, Spectra, and Redox Properties of Myeloperoxidase -A Density Functional Study[END_REF][START_REF] Battistuzzi | Influence of the Covalent Heme-Protein Bonds on the Redox Thermodynamics of Human Myeloperoxidase[END_REF] . The high E °' values were suggested to form the thermodynamic basis that enables Cl ➖ oxidation by MPO [START_REF] Furtmuller | Standard Reduction Potentials of All Couples of the Peroxidase Cycle of Lactoperoxidase[END_REF][START_REF] Arnhold | Redox Properties of the Couple Compound I/Native Enzyme of Myeloperoxidase and Eosinophil Peroxidase[END_REF] .

Table 2 Standard reduction potentials (pH 7) of all relevant redox couples of MPO, mutant MPO, LPO and horseradish peroxidase (HRP) as example of heme b -containing peroxidases having a histidine as proximal ligand, but not featuring the covalent posttranslational modifications of MPO and LPO. (This table is constructed from Table 3 of the review of Nicolussi et al. [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF] . We also point to the references given heir in) Based on crystallographic data, it has been argued that the covalent bonds formed between the tetrapyrrole ring and the protein, especially the sulfonium linkage, induce significant distortions of the heme from planarity, [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF][START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF][START_REF] Fenna | Structure of the Green Heme in Myeloperoxidase[END_REF] and that this plays a role in determining the reactivity. The first crystal structure of MPO obtained to 1.8 Å resolution has led to the observation that in MPO, the heme group features a pronounced bow-shape structure with low symmetry [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF][START_REF] Fiedler | X-Ray Crystal Structure and Characterization of Halide-Binding Sites of Human Myeloperoxidase at 1.8 A Resolution[END_REF][START_REF] Fenna | Structure of the Green Heme in Myeloperoxidase[END_REF] . This pronounced out-of-plane distortion has been suggested to cause the different spectral and redox properties of MPO compared to LPO [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF] . It is well established that hemes in proteins rarely adopt planar conformations, which are generally conserved within functional classes of proteins [START_REF] Shelnutt | Nonplanar Porphyrins and Their Significance in Proteins[END_REF] and that modulation of ruffling regulates electron transfer with implications for biological activity [START_REF] Sun | Investigations of Heme Distortion, Low-Frequency Vibrational Excitations, and Electron Transfer in Cytochrome c[END_REF] . For an overview, see Kleingardner et al. [START_REF] Kleingardner | Biological Significance and Applications of Heme c Proteins and Peptides[END_REF] However, the suggestion that the distortion of the prosthetic group by the sulfonium linkage in MPO gives rise to the highly positive E °' values is in contradiction with other studies based on DFT calculations and NMR experiments on ferric Fe(+III) heme. These studies conclude the decrease of the heme reduction potential in response to increased ruffling [START_REF] Liptak | NMR and DFT Investigation of Heme Ruffling: Functional Implications for Cytochrome c[END_REF][START_REF] Kleingardner | The Influence of Heme Ruffling on Spin Densities in Ferricytochromes c Probed by Heme Core 13 C NMR[END_REF] . Out of plane deformations introduce several changes in the electronic structure of the low-spin His/Met ligated ferric heme cofactor: firstly, the positive spin densities on the pyrrole carbon decrease as a function of ruffling. Secondly, the negative spin density on the heme meso carbons decreases in ruffled hemes. Thirdly, all the occupied Fe 3d-based molecular orbitals are destabilized by ruffling, with the Fe-3dxy based MO being the most sensitive. These changes are consistent with the lowering of the Fe(+III)/Fe(+II) reduction potential.

E °' (
These ambiguities and fragmentary studies present in the current literature towards the oxidation of Cl ➖ by MPO encouraged us to re-open the investigation to the oxidizing power of MPO.

Computational details

In total, nine models (see Scheme 1, giving an overview) have been constructed to mimic HRP, MPO, LPO and CPO starting from available high-resolution crystal structures with PDB code 5fiw (for the MPO and LPO models) and 2j5m (for the CPO model). For each model, the geometry of all possible spin states of all redox intermediates (resting state, Cpd 0, Cpd I, Cpd II and protonated Cpd II) and of the reactant state (Cpd I.H2O.Cl ➖ ), the transition state (TS) and the product complex (PC) for the oxidation of Cl ➖ have been fully optimized. IRC calculations show a continuous path from RS, over TS to PC. All geometries were obtained with Gaussian09 35 at the B3LYP/6-31G(d) level in an SMD solvent model using diethyl ether (ε=4.2400) as solvent. For the MPO imi and HRP imi models, the geometries of the doublet RS, TS and of Cpd II, protonated Cpd II and of the resting state are obtained with Gaussian16 36 . Several tests confirmed that both versions of the calculation software give similar geometries (atomic distances are the same up to 0.01 Å), varying only a few hundreds to a (few) tenth of one kcal/mol in energy. Additional single point calculations have been performed at the B3LYP/6-311+(2d,p) level in an SMD solvent model using diethyl ether (ε=4.2400) as solvent, with Gaussian16 36 . The empirical dispersion correction for density functionals (D3) with Becke and Johnson (BJ) damping has been added during geometry optimization and during single point calculations [START_REF] Grimme | A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu[END_REF][START_REF] Grimme | Effect of the Damping Function in Dispersion Corrected Density Functional Theory[END_REF] . B3LYP is commonly used and generally accepted as reliable functional [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF] for heme systems and found to fairly well reproduce multireference calculations [START_REF] Altun | Multireference Ab Initio Quantum Mechanics/Molecular Mechanics Study on Intermediates in the Catalytic Cycle of Cytochrome P450cam[END_REF][START_REF] Radoń | Spin-State Energetics of Heme-Related Models from DFT and Coupled Cluster Calculations[END_REF] . Free energies G 0 were obtained as the sum of the electronic energy (Eelec), the zero point energy (ZPE) and a series of thermal correction terms (G 0 = Eelec + ZPE+ Erot + Evib + Etrans + RT -TS; with S: entropy, R: gas constant, T: temperature and Erot , Evib, Etrans: the rotational, vibrational and translational energy), calculated at a standard state of 1 M H + and 298.15 K. Eelec was calculated at the B3LYP/6-311+(2d,p) D3 BJ level in an SMD solvent model using diethyl ether (ε=4.2400) as solvent. ZPE and thermal correction terms were calculated from partition functions using the harmonic oscillator approximation, in a frequency calculation following the geometry optimization, at the same calculation level. To correct for the breakdown of the rigid rotor harmonic oscillator (RRHO) approximation as approach to calculate the vibrational entropy, all vibrational frequencies below the cut-off value of 50 cm -1 were uniformly shifted to 50 cm -1 (quasi harmonic approximation as suggested by Cramer and Truhlar et al. [START_REF] Ribeiro | Use of Solution-Phase Vibrational Frequencies in Continuum Models for the Free Energy of Solvation[END_REF] ). Reaction free energies ' G  at pH=7 are obtained from the calculated G 0 values plus 7 ( ' ') ln10 p r RT --, making the conversion between the standard state of 1 M H + (pH=0) and pH=7. ' p and ' r are the amounts of H + ions respectively formed and used during the reaction.

Results and discussion

To understand the chemistry of Cl ➖ oxidation by Cpd I for the various heme enzyme models (cfr. eq. ( 4)), we have studied the transition state, some of the intermediates, and their relevant thermodynamic properties in order to get an overview of the reaction mechanisms, and the factors that determine the selectivity between them. We have chosen a broad range of models (Scheme 1) approximating the heme containing active site of HRP, MPO, LPO and CPO. By choosing a cluster model and the semi-quantitatively accurate B3LYP DFT approach, our study cannot yield quantitative results, but has been designed to provide insights into the chemistry of these peroxidases.

- - HRP Imidazole -CH3 -CH3 -C2H3 HRP imi Imidazolate -CH3 -CH3 -C2H3 MPO-E242Q Imidazole -CH3 -CH3 -C2H2Met243 LPO Imidazole -CH2Asp94 -CH2Glu242 -C2H3 MPO Imidazole -CH2Asp94 -CH2Glu242 -C2H2Met243 MPO imi Imidazolate -CH2Asp94 -CH2Glu242 -C2H2243 CPO -CH3S ➖ Cys29 -CH3 -CH3 -C2H3 CPO envir -CH3S ➖ Cys29, Ala31, Leu32 PDB numbering 2j5m -CH3 -CH3 -C2H3
In the minimal model, the porphyrin ring of the heme group is used without any side-chains, i.e. with with hydrogen atoms replacing the propionate and vinyl groups and/or their covalentlymodified forms present in some of the peroxidases. The proximal histidine ligand is replaced by a (neutral) imidazole group. This core model serves as a basis for all the other models. Its complexity is gradually extended, first by adding the vinyl-groups and -CH3 groups replacing the propionate groups (HRP model). Secondly, the two ester linkages that are present in LPO and MPO were added (LPO model) and thirdly, the (positively charged) sulfonium linkage is added (MPO model). For CPO, an HRP-like model is constructed, but with the proximal cysteine replaced by (a negatively charged) thiolate.

To investigate the suggested deprotonation of the proximal histidine [START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF] , for both the HRP and MPO models, a second model was built with a negatively charged imidazolate as proximal ligand, the HRP imi and MPO imi models.

Coordinates of the optimized geometries of the redox states participating in the reaction pathways can be found in SI, Table S1. The energetic order of the spin states of Cpd 0, Cpd I and of the resting state has been validated against available literature (Table S2) [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF][START_REF] Rydberg | On the Role of the Axial Ligand in Heme Proteins: A Theoretical Study[END_REF] . The lowest energy spin state of each redox intermediate is used to obtain the thermodynamic properties and activation energies, which is the doublet state for Cpd 0 and the quartet/doublet state for the Fe(III) resting state when respectively imidazole or thiolate is present as proximal ligand.

Cpd II and its protonated form adopt a triplet ground spin state [START_REF] Lai | What Kinds of Ferryl Species Exist for Compound II of Chloroperoxidase? A Dialog of Theory with Experiment[END_REF][START_REF] Derat | Two-State Reactivity, Electromerism, Tautomerism, and "Surprise" Isomers in the Formation of Compound II of the Enzyme Horseradish Peroxidase from the Principal Species, Compound I[END_REF] . For the degenerate doublet/quartet state of Cpd I, the quartet state has been chosen (though some additional reactivity calculations were performed for the doublet state).

Reactant state formation

Structure Preceding geometry optimization of the reactant complex (Cpd I.H2O.Cl ➖ ), H2O.Cl ➖ was manually positioned at the distal side of the heme. In the optimized geometry, Cl -accepts a hydrogen bond from the water molecule spanning Cl ➖ and the oxygen atom of Cpd I. Cl ➖ is positioned above a meso carbon, at a distance of 2.8 Å to the meso carbon, at 4.8 Å to Fe and at 4.6 Å to the oxygen of Cpd I (distances reported for the optimized MPO model).

In structure 1d2v [START_REF] Fiedler | X-Ray Crystal Structure and Characterization of Halide-Binding Sites of Human Myeloperoxidase at 1.8 A Resolution[END_REF] , showing bromide binding to the resting state of MPO in a distal cavity 3.5-3.6 Å away from N of the distal histidine, Br ➖ takes a similar position, with a Br ➖ ---meso carbon distance of 5.6 Å and a Br ➖ ---Fe (+III) distance of 5.0 Å. All other halides present in structure 1d2v and the halides present in all other available PDB structures are found in cavities remote from the heme [START_REF] Fiedler | X-Ray Crystal Structure and Characterization of Halide-Binding Sites of Human Myeloperoxidase at 1.8 A Resolution[END_REF] . Also for CPO, it has been suggested, based on measured ferric-ferrous redox potentials that are not affected by halide binding [START_REF] Makino | Oxidation-Reduction Potential Measurements on Chloroperoxidase and Its Complexes[END_REF] , that halides do not bind as axial ligand to the ferric heme. As such, halides might first bind into remote binding sites before going to the heme pocket. This questions the assumed reactant complex, especially for CPO, and might lead to erroneous (high) H2O.Cl ➖ -binding energies calculated in our model system approach (Table 3).

In the given models, a small favorable effect of the ester linkages on Cl ➖ binding is observed.

No direct experimental information on Cl ➖ binding to Cpd I is available, but the calculated effect of the ester groups agrees with experimental data on cyanide binding to the MPO resting state. The dissociation constants KD obtained for cyanide binding from respectively kinetic experiments and by spectroscopic titration are 1.9 M and 1.4 M for wild type recMPO and 2.4 M and 3.5 M for the Glu242Gln mutant, measured at pH=7 and 25°C [START_REF] Zederbauer | Role of the Covalent Glutamic Acid 242-Heme Linkage in the Formation and Reactivity of Redox Intermediates of Human Myeloperoxidase[END_REF] .

H2O.Cl -binding seems to be strongly favored in the presence of the sulfonium linkage and strongly disfavored by the presence of the negatively charged thiolate and imidazolate proximal ligands (Table 3). This strong unfavorable binding of Cl -to CPO is in line with the experimental finding that halides do not bind on axial positions [START_REF] Makino | Oxidation-Reduction Potential Measurements on Chloroperoxidase and Its Complexes[END_REF] . Some caution is however necessary as the calculated strong effects most probably arise from the approximated models not correctly describing the electrostatic interactions around the ligand's binding.

Electronic structure Cpd I is found to have a 'classic' Cpd I electron configuration, in the models for CPO, MPO and MPO imi that we considered. To show this, we considered orbital coefficients for half-filled orbitals and the Mulliken spin densities (Table 4). The 'classic' electronic structure of Cpd I also extends to the situation where H2O.Cl -is bound in the reactant complex. This standard configuration for Cpd I involves the presence of three unpaired electrons, of which two are situated on the Fe-O center, residing in two singly occupied p * FeO , p * xz and p * yz , orbitals, involving anti-bonding interactions between the metal 3d xz , yz and the oxygen 2 p x,y orbitals [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF] . The third unpaired electron is located in a high-lying non-bonding porphyrin orbital. The porphyrin's singly occupied orbital of Cpd I do not mixes with the imidazole/imidazolate orbitals and has merely an A 1u character (Figure 3), similar to what was described before for horseradish peroxidase (HRP) [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF] . Due to the small mixing of the imidazole/imidazolate orbitals with a 2u , the gap between a 2u and a 1u remains small and thus the gap between the A 1u and A 2u states remains small. The heme porphyrin in the optimized models with imidazole/imidazolate as proximal ligand is fairly largely distorted from planarity; a distortion that is remarkably more compared to what we observe from our optimized geometries when the charged thiolate is present. These out-of-plane distortions lowers the symmetry of the porphyrin so that the A 1u and A 2u states can mix their a 2u and a 1u orbitals and hence also their characters.

Upon Cl -binding to Cpd I, the single occupied porphyrin orbital becomes more localized and is as such a more pronounced a 2u orbital. This can also be seen from the localization of the Mulliken spin densities on the meso carbon atoms (Table S3), despite a similar out-of-plane distortion of the porphyrin ring of the Cpd I.H2O.Cl -complex compared to Cpd I. The positive spin density present on the substituent vinyl groups of the porphyrin is also removed upon Cl - binding (Table S3). It is known that many close-lying valence states with different charge distributions, contribute to the electronic structure of Cpd I species [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF][START_REF] Derat | Principal Active Species of Horseradish Peroxidase, Compound I: A Hybrid Quantum Mechanical/Molecular Mechanical Study[END_REF][START_REF] Lai | Effects of Substrate, Protein Environment, and Proximal Ligand Mutation on Compound I and Compound 0 of Chloroperoxidase[END_REF] . We assume that the chloride by interacting differently with these valence states perturbs the electronic structure.

For CPO, the singly occupied porphyrin orbital strongly mixes with the donating thiolate ligand in CPOas was described before for P450 [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF] . This mixing is less present when Cl -is bound.

For the CPO model, a drop of the Mulliken spindensity on the CH3S -ligand from 0.695 to 0.325 is seen, while the spindensity on the porphyrin ring increases from 0.201 to 0.761 upon H2O.Cl - binding.

Transition state and reaction path

As well as exploring the electronic structure of Cpd 0, I and II, and the thermochemistry of oneand two-electron reaction steps (see next paragraph), we also wish to carry out a preliminary study of the reactivity towards halide (and in particular chloride) oxidation. To do this, transition states (TS) for Cl -oxidation are obtained for the different heme models and reaction free energies are reported in Table 3. In principle, such oxidation TSs could exist on both the quartet and the doublet potential energy surface. Previous work on analogous two-electron oxidation of sulfide [START_REF] Rydberg | Sulfur, and Nitrogen Oxidation and Dealkylation by Cytochrome P450[END_REF] and amines [START_REF] Lonsdale | Quantum Mechanics/Molecular Mechanics Modeling of Drug Metabolism: Mexiletine N-Hydroxylation by Cytochrome P450 1A2[END_REF] suggests that the doublet TS will lie lower in energy (for oxidation of alkyl halides, however, the quartet TS seems to be lower [START_REF] Ji | Metabolism of Halogenated Alkanes by Cytochrome P450 Enzymes. Aerobic Oxidation versus Anaerobic Reduction[END_REF] ). In this preliminary study, however, although we have located the doublet TS for one model, we have focused mostly on the quartet TSs, which already provide the required insight into reactivity trends. The doublet TSs proved harder to locate, perhaps because they seem to be low in energy and due to shortcomings in the treatment of the solvation. As previously mentioned, our cluster model and DFT approach precludes reaching quantitative results and seems to affect our ability to study the doublet TS. This also means that our results for reactivity are not quantitative, although trends should be more reliable. For the quartet TSs, a similar structure and electronic structure pattern is found for all of the models (Figure 4). As such, neither the charge of the proximal ligand nor its nature (imidazole or thiol) nor the heme substituents (covalent ester and sulfonium linkages) seem to influence the nature of the transition state. S4) and orbital coefficients point to a productlike TS. The positive spin density on the porphyrin ring and on the oxygen atom, present in the Cpd I reactant state, is removed, while the spin density on Fe increases (compare Table 4 andS4). The doublet MPO imi TS is the lowest energy state, compared to the quartet state. Also, for the product complex, a doublet spin state as ground state is found. As such, the degeneracy between the doublet and quartet spin state found in the reactant state is lifted when evolving to the product complex through the TS. For the oxidation of amines to hydroxylamines [START_REF] Lonsdale | Quantum Mechanics/Molecular Mechanics Modeling of Drug Metabolism: Mexiletine N-Hydroxylation by Cytochrome P450 1A2[END_REF] , the doublet reaction path is also the lowest energy-path, with the degenerate reactant state evolving to a doublet product state through a doublet TS. For the C-H hydroxylation, the doublet and quartet spin state profiles are very close in energy [START_REF] Shaik | Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes[END_REF] . An IRC (Figure S1a) confirms the continuous path from reactant state to transition state. The analysis of the Mulliken charges (Figure S1b) along the IRC confirms some charge transfer from Cl -before reaching the TS, although the majority of the charge is transferred after the TS.

Free energy of activation Table 3 reports the calculated free energies of activation G ‡ for the Cl ➖ oxidation by Cpd I in the given models, which can be related to experimental k cat values. These are not directly available from experiment for the key halide oxidation step, but can be approximately extracted from analysis of the apparent second order rate constants are (Table 1). From these data, we derived the corresponding k cat values (see SI Text S1 for the details of the derivation). It should also be noted that the comparison with experiment is difficult to make in a quantitative way due to the fact that our models do not provide a definitive treatment of the doublet vs. quartet TSs, and the fact that for some of the models, very strong electrostatic interactions between the chloride and the Cpd I model influence the calculated energies significantly. A rough description of these interactions is included through the use of a microsolvating water molecule for the chloride reactant, and through use of a continuum solvent model. However, based on the large magnitude of the differences in calculated free energies for heme models with different total charge, it is likely that a more accurate treatment would be required in order to obtain even semi-quantitatively accurate free energies. Such a treatment goes beyond the scope of this work, which examines instead the comparative properties of the different peroxidase models. It is not clear if a reactant state complex is formed and thus both the reactant complex Cpd I.H2O.Cl ➖ or the unbound Cpd I + H2O.Cl -could be valid choices to take as reference state to calculate the reaction barriers. Reaction barriers relative to both reference states are reported in Table 3 and in the following sections, the barriers with the reactant complex as reference state are discussed. A quantitative discussion of the reactivity would require a more careful investigation of this point.

For CPO, apparent second order rate constants are obtained at pH=4.8 (Table 1). The binding of nitrate, acetate and formate is documented at pH 3 and the dissociation constants are measured to be respectively 0.55 mM, 93 mM and 0.135 mM [START_REF] Sono | Ligand and Halide Binding Properties of Chloroperoxidase: Peroxidase-Type Active Site Heme Environment with Cytochrome P-450 Type Endogenous Axial Ligand and Spectroscopic Properties[END_REF] . Titration results show that halides bind in their anionic form to the acid form of the ferric (Fe +III) enzyme [START_REF] Sono | Ligand and Halide Binding Properties of Chloroperoxidase: Peroxidase-Type Active Site Heme Environment with Cytochrome P-450 Type Endogenous Axial Ligand and Spectroscopic Properties[END_REF] . Based on the size of chloride compared to nitrate/formate versus acetate, we assume a KM or Cl -of 0.1 mM. Given the experimentally measured [START_REF] Zhang | Kinetics of Two-Electron Oxidations by the Compound I Derivative of Chloroperoxidase, a Model for Cytochrome P450 Oxidants[END_REF]

k app = k cat K M of 2.5 10 4 M -1 s -1 , a
k cat value of 2.5 s -1 can be assumed. This corresponds to an estimated experimental free energy of activation G ‡ of 16.9 kcal/mol. For CPO, the existence of a reactant complex is questionable and thus we end up with a calculated activation energy of G ‡,4 = 44.52 kcal/mol, obviously too high. Optimizing a TS in a model that includes the hydrogen bonding environment around the proximal thiolate was not possible, but the energy values point to a decrease of the reaction barrier when the negative charge is shielded by the environment (Table 3b).

As the binding of chloride to spleen MPO has been found to be pH dependent [START_REF] Ikeda-Saito | A Study of Ligand Binding to Spleen Myeloperoxidase[END_REF] , the experimentally obtained apparent second order rate constant k app' equals :

k app' = k cat K M K A H 3 O + é ë ù û (5)
measured at pH = 7.0 and 5.0 for the Cl ➖ oxidation by MPO (Table 1). Data for KM are very rare, but an estimate of KM for Cl -binding to Cpd I can be made from the available data of Cl - binding to the resting state of MPO. At room temperature (20° C) [START_REF] Ikeda-Saito | A Study of Ligand Binding to Spleen Myeloperoxidase[END_REF] , KCl -is reported to be 0.06 mM; with a pKa of the resting state of 3.9. Assuming that binding to the resting state and Cpd I gives approximately the same values, we can approximate the dissociation constant KM for chloride as 0.05 mM and KA as 10 -4 M. With [H3O + ] at pH 7 = 10 -7 M and k app' of 2.5 10 4 M -1 s - 1 (Table 1), k cat equals 0.125 10 4 s -1 . For the E242Q mutant, KM and KA are assumed to have the same value as for wild type MPO and thus with k app' of 0.0065 10 4 M -1 s -1 (Table 1), an estimated k cat of 0.0003 10 4 s -1 at pH= 7.0 is obtained. These estimated k cat values can now be compared directly to the calculated activation free energies G ‡ , using the reversed Eyring equation. For MPO, the experimentally estimated G ‡ equals 13.2 kcal/mol; for the E242Q mutant, G ‡ equals 16.8 kcal/mol, to which the calculated G ‡ values (Table 3) largely agrees: for MPO, G ‡,4 = 11.76 kcal/mol and for MPO-E242Q G ‡,4 = 14.95 kcal/mol. Both values are lower than the G ‡,4 of 16.08 kcal/mol calculated for the HRP model, consistent with the experimentally measured oxidation activity (Table 1). The calculated G ‡,4 values go together with the electron withdrawing nature of the ester and sulfonium, hereby enhancing the oxidizing power of Cpd I. This is consistent with the increasing/decreasing Mulliken spin density on Cpd I's oxygen/Fe in the presence of the esters and sulfonium (Table 4), allowing the substituted heme to easily accept electrons. Switching to a charged imidazolate as proximal ligand (HRP imi and MPO imi models), raises the reaction barrier in the quartet path slightly, to 17.73 kcal/mol for the HRP model, and to 14.40 kcal/mol for the MPO model (Table 3). The doublet path for the MPO imi model lowers the barrier to 4.65 kcal/mol. As for CPO, the existence of the assumed reactant state for MPO imi is also questionable. With the unbound Cpd I + H2O.Cl ➖ state as reference, the reaction barriers of the doublet and quartet path calculated for the MPO imi model are equal (11 kcal/mol), in agreement with the estimated experimental value of 13.2 kcal/mol.

Thermodynamic properties of reaction intermediates involved in Cl ➖ oxidation

As depicted in Figure 1, Cpd I can be reduced back to the native ferric state following two competing paths, a one and a two-electron path. A rough estimate of the reactivity of Cpd I in the one-and two-electron cycles comes from the calculation of thermodynamic properties as the one-and two-electron affinity of Cpd I and the proton affinity of Cpd II (Table 5). These calculated properties are approximations for the experimentally measured one-and two-Cpd I reduction potentials (Table 2). As described in the previous paragraph, approximated models are not likely to describe well the electrostatics around the active site and thus, some caution is necessary when discussing the electron affinities. The very negative one-electron affinities of Cpd I in the presence of the sulfonium linkage (MPO E242Q and MPO model) (Table 5) is however consistent with the highly positive E °' value of the Cpd I/Cpd II redox couple, experimentally measured to be 1350 mV [START_REF] Furtmüller | Redox Properties of the Couples Compound I/Compound II and Compound II/Native Enzyme of Human Myeloperoxidase[END_REF] . The ester groups have only a minor influence on the one-electron affinity, which can't directly be linked to the experimental data of Table 2. The calculated one-electron affinity of Cpd I in the CPO and in the imidazolate bearing models is less negative, indicating a lower one-electron affinity. Also, the presence of the negatively charged chloride ligand, diminishes the one-electron affinity.

These calculated electron affinities follow the sequence of the stability of Cpd I, calculated as the energy difference between Cpd I and Cpd 0 [START_REF] Lai | Effects of Substrate, Protein Environment, and Proximal Ligand Mutation on Compound I and Compound 0 of Chloroperoxidase[END_REF] . The two electron reduction potential is much less affected by the sulfonium linkage as is suggested from calculated redox potentials [START_REF] Devarajan | Effect of Covalent Links on the Structure, Spectra, and Redox Properties of Myeloperoxidase -A Density Functional Study[END_REF] and from the available E °' values of the Cpd I/Fe(III) redox couple (Table 2) In agreement with these data, our calculated two-electron affinities of Cpd I in the different models vary less by the presence of the sulfonium linkage than the one-electron affinities (Table 5). The two-electron affinities are calculated to be positive and thus substantially more unfavorable that the one-electron affinities.

One important feature, extensively discussed [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF][START_REF] Yosca | Iron(IV)Hydroxide PK(a) and the Role of Thiolate Ligation in C-H Bond Activation by Cytochrome P450[END_REF][START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] as main mechanism to suppress one-electron oxidations by Cpd I, is the pKa of Cpd II. The electron donating proximal thiolate and tyrosinate ligands elevate the pKa of Cpd II as compared to neutral ligands such as histidine. This effect is often referred to as the push-effect. Cpd II's pKa can be estimated by the proton affinity (PA) of Cpd II calculated as the energy difference between its protonated and un-protonated form. The calculated Cpd II PA in models having a neutral histidine as ligand has a substantially less negative value and is thus lower than the Cpd II PA calculated in models having a charged thiolate as ligand. This agrees with experimentally measured Cpd II pKa values of £ 3.5 in proteins having a neutral histidine proximal ligand [START_REF] Yosca | Setting an Upper Limit on the Myoglobin Iron(IV)Hydroxide PKa: Insight into Axial Ligand Tuning in Heme Protein Catalysis[END_REF][START_REF] Behan | On the Status of Ferryl Protonation[END_REF] , while for thiolate and tyrosinated ligated proteins, a high Cpd II pKa of about 12 (P450) -13 (catalase) has been reported [START_REF] Green | Bond Activation in Heme Proteins: The Role of Thiolate Ligation in Cytochrome P450[END_REF][START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] . Installing an imidazolate as proximal ligand shifts the Cpd II PA towards more negative values, resembling as such more the Cpd II PA calculated for the CPO models.

Table 3: Reaction free energies

' G   (kcal/mol) (a) and reaction energies (b) calculated for the different models with proximal imidazole, imidazolate and thiolate, for the quartet reaction path. 𝐺° of the small molecules present in the reaction equations are documented in Table S5. 4) on the oxygen/Fe atoms of its Cpd I, allowing the substituted heme to easily accept electrons, and thus in a low reaction barrier for Cl ➖ oxidation (Table 3). As depicted in Figure 1, Cpd I can be reduced back to the native ferric state following two competing paths, a one-and a two-electron path. Having a neutral proximal histidine, MPO's Cpd I one-electron affinity is very high (strong negative value, relative to the calculated values for the other models), which is accompanied with a low Cpd II's proton affinity (small negative value, relative to the calculated values for the other models) (Table 5).

a ' G   (kcal/mol) Minimal model HRP HRP imi MPO - E242Q LPO MPO MPO imi MPO imi CPO $=4 $=2 $=4
As such, compared to CPO, MPO's one-electron pathway (Figure 1) seems to be strongly favored over the two-electron oxidation cycle. The question remains how MPO prevents going into these futile one-electron cycles and is as such able to efficiently catalyze the two-electron oxidation of Cl ➖ . Surprisingly, the answer might lie in the installation of an imidazolate as proximal ligand.

The most obvious choice for the protonation state of the proximal histidine ligand in MPO seems to be the neutral imidazole form. This reasoning is based on earlier work (including QM/MM calculations, Raman spectroscopy and mutagenesis studies) on cytochrome c peroxidase (CcP) [START_REF] Smulevich | Heme-Protein Interactions in Cytochrome c Peroxidase Revealed by Site-Directed Mutagenesis and Resonance Raman Spectra of Isotopically Labeled Hemes[END_REF][START_REF] Smulevich | Fifteen Years of Raman Spectroscopy of Engineered Heme Containing Peroxidases: What Have We Learned?[END_REF] , ascorbate peroxidase (APX) [START_REF] Bathelt | QM/MM Studies of the Electronic Structure of the Compound I Intermediate in Cytochrome c Peroxidase and Ascorbate Peroxidase[END_REF][START_REF] Harvey | QM/MM Modeling of Compound I Active Species in Cytochrome P450, Cytochrome C Peroxidase, and Ascorbate Peroxidase[END_REF] and HRP [START_REF] Derat | Principal Active Species of Horseradish Peroxidase, Compound I: A Hybrid Quantum Mechanical/Molecular Mechanical Study[END_REF] , indicating the imidazole rather than the imidazolate nature of the proximal His. The former proteins have a proximal histidine ligand forming a strong hydrogen bond connecting the imidazole's -NH group with the oxygen of an aspartate. In MPO, LPO and EPO, the environment around the proximal His is different. Here, the proximal His is involved in a hydrogen-bond network (Figure 5) with an Asn and Arg, all at Arg-NH2---OAsn and Asn-NH2---NHis distances < 3.1 Å, comparable to what has been reported by Carpena, X. et al [START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF] . Structural data and molecular modelling show the -NH2 group of Asn interacting with N  of the proximal His, whereas its carbonyl group interacts with Arg at short H-bonding distances. It has been argued that a neutral proximal His would install unfavorable contacts [START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF][START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF] at these short distances. In MD simulations with imidazole, the His-Asn-Arg H-bonds are substantially longer. Short His-Asn-Arg hydrogen bond distances are found in all checked crystal structures of MPO and LPO available in the PDB. For EPO, no experimental structures are known, but using homology modelling (for details, see SI text S2), we arrived at a structural arrangement around the proximal His similar to this of MPO and LPO (Figure 5). Also, the FeHis stretching mode found in MPO and LPO at fairly high frequencies (248-244 cm -1 ) might be consistent with a substantial imidazolate character of the proximal ligand in these proteins.

The chemistry of the oxidation reaction is seemingly not influenced by the imidazole/imidazolate character of the proximal ligand: a 'classic' Cpd I electron configuration is found and the (electronic) structures of the TS are similar to these obtained for the other models (Figure 3). Ligation with an electron donating anion -either thiolate, as for example in cytochrome P450 [START_REF] Yosca | Iron(IV)Hydroxide PK(a) and the Role of Thiolate Ligation in C-H Bond Activation by Cytochrome P450[END_REF] and chloroperoxidases, or tyrosinate as in catalases [START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] has been pointed forward to be crucial for Cpd I to effectively catalyze two-electron oxidations. The electron donating proximal thiolate and tyrosinate ligands elevate the pKa of Cpd II as compared to neutral ligands as histidine. A basic iron(IV)oxo species, installed by anion ligation, minimizes off-pathway chemistry by elevating the pKa of Cpd II, hereby lowering the Cpd I one-electron reduction potential [START_REF] Yosca | Spectroscopic Investigations of Catalase Compound II: Characterization of an Iron(IV) Hydroxide Intermediate in a Non-Thiolate-Ligated Heme Enzyme[END_REF] . The relative free energy towards the productive two-electron pathway is as such shifted. Based on Marcus theory, this will suppress the rate constant for one-electron oxidations and thus futile one-electron cycles (Figure 1) are avoided [START_REF] Yosca | A New Look at the Role of Thiolate Ligation in Cytochrome P450[END_REF][START_REF] Yosca | Iron(IV)Hydroxide PK(a) and the Role of Thiolate Ligation in C-H Bond Activation by Cytochrome P450[END_REF] . Within the limitations of approximate models, the calculated one-electron affinity of Cpd I and proton affinity of Cpd II (Table 5) seems to point to a similar mechanism for CPO; i. e. disfavoring one-electron chemistry by elevating the proton affinity of Cpd II (as indicator for the Cpd II pKa) and lowering the Cpd I electron affinity (as indicator for the one electron reduction potential of Cpd I). These favorable thermodynamics are also found for the imidazolate bearing models as compared to models with imidazole (Table 5). This can strongly suggest the imidazolate nature of the proximal ligand in MPO, but is on the other hand contradictory to the high one-electron reduction potentials measured for the Cpd I/Fe(III) and the Cpd I/Cpd II redox couples (Table 2). The reaction barriers (when the reactant state complex is taken as reference state) are not influenced strongly by the imidazolate/imidazole nature of the proximal ligand. When the unbound Cpd I + H2O.Cl ➖ state is taken as refence state, the reaction barriers with imidazolate are substantially higher than with imidazole as proximal ligand. This reflects most probably the poor description of the Cl ➖ binding with the electrostatic effect being too large in the imidazolate models. Fully solvated models will be needed to describe correctly the electrostatics. Further, inclusion of specific solvent molecules in the vicinity of the heme will shed light on their possible roles in proton transfer pathways [START_REF] Kwon | Visualizing the Protons in a Metalloenzyme Electron Proton Transfer Pathway[END_REF] .

The calculated G ‡ values for all models having an imidazole or imidazolate as proximal ligand are fairly low. Especially the sulfonium linkage has an activating effect, while the ester linkages play only a small activating role. Therefore, LPO and EPO not bearing the sulfonium linkage, but only the ester linkages cannot oxidize Cl ➖ (LPO) or only to some extent (EPO), but still can oxidize Br ➖ , CN ➖ and I ➖ , having a substantially lower standard reduction potential. HRP, not having the ester linkages does not exert halide oxidation. The calculated G ‡ values (Table 3) agrees with this sequence of activity of the oxidation reaction.

Our study, consistent with previous work by Ryde et al [START_REF] Devarajan | Effect of Covalent Links on the Structure, Spectra, and Redox Properties of Myeloperoxidase -A Density Functional Study[END_REF] , could not find evidence that out-ofplane distortion of the heme, as has been previously suggested [START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF][START_REF] Carpena | Essential Role of Proximal Histidine-Asparagine Interaction in Mammalian Peroxidases[END_REF][START_REF] Fenna | Structure of the Green Heme in Myeloperoxidase[END_REF] , is at the basis of the ability of MPO to oxidize Cl ➖ . No convincingly more out-of-plane distortion of the heme porphyrin in the presence of the sulfonium linkage could be observed from visual inspection of the dihedral angles in the porphyrin ring of the optimized geometries of the different model systems.

Obtaining a clear view in the approximated models on the role of Cl ➖ binding in optimizing the catalytic efficiency (

k cat K M
) was also not possible using approximate models. A model for substrate binding to Cpd I has been made based the binding of halides to the cyanide complex of MPO. Energetically similar positions for H2O.Cl ➖ at the distal side of the heme pocket can be chosen as starting position, giving raise to structural uncertainties. Also, electrostatic effects might not be fully covered when using approximate models. Electron paramagnetic resonance (EPR) spectra suggest the disruption of the original distal hydrogen bond network in the MPO-CN complex and the involvement of the distal histidine and arginine residues. These data also suggest direct contact between the halide and the iron sixth ligand at low pH only [START_REF] Ikeda-Saito | A Study of Ligand Binding to Spleen Myeloperoxidase[END_REF][START_REF] Lee | Interaction of Halides with the Cyanide Complex of Myeloperoxidase: A Model for Substrate Binding to Compound I[END_REF] .

The existence of the Cpd I.H2O.Cl ➖ reactant complex might be questionable. This is especially true for CPO for which it has been experimentally suggested that halides do not bind on axial positions [START_REF] Makino | Oxidation-Reduction Potential Measurements on Chloroperoxidase and Its Complexes[END_REF] . The experimentally obtained KM data might thus refer to binding in remote pockets and not to heme binding directly. We anticipate this to be a second order effect with an overall impact of going from a remote binding site to the heme pocket on the activation energy presumably being small and based on available experimental data [START_REF] Ikeda-Saito | A Study of Ligand Binding to Spleen Myeloperoxidase[END_REF][START_REF] Sono | Ligand and Halide Binding Properties of Chloroperoxidase: Peroxidase-Type Active Site Heme Environment with Cytochrome P-450 Type Endogenous Axial Ligand and Spectroscopic Properties[END_REF] , we have assumed a KCl - in the mili-molar range for both CPO and MPO to estimate the kcat, based on experimentally measured apparent second order rate constants (Table 1) .

To fully judge the effect of out-of-plane distortions and to obtain a clear view on Cl ➖ binding, studies on entire protein systems should be undertaken.

Conclusion

Although further studies on entire protein systems will be useful to correctly describe the electrostatics and to fully grasp the contribution of ligand binding and heme dynamics on reactivity, we can already conclude from this model-system based study that the histidine bearing MPO has everything to be an as good oxidant as CPO, having a proximal cysteine thiolate. The ester and especially the sulfonium linkages are hiking up the reactivity of MPO to enable chloride oxidation. The imidazole/imidazolate character of the proximal ligand does not influence the chemistry of the oxidation reaction, but an imidazolate proximal ligand, installed by the H-bonding environment involving a conserved Asp and Arg, instead of a neutral His seems to be necessary to shift the thermodynamic equilibrium from a sequential one-electron to a two-electron event.

i+ 3 )Figure 1

 31 Figure 1 Schematic representation of the one-and two-electron pathways of MPO Reaction 1: ferric MPO (Fe +III; resting state) is oxidized by hydrogen peroxide to compound I (i.e., oxoiron(IV) porphyrin radical cation). Reaction 2: compound I is reduced to compound II (i.e., protonated oxoiron(IV)) by a one-electron donor. Reaction 3: compound II is reduced to the native state, thereby oxidizing a second substrate molecule. Reactions 1, 2, and 3 constitute the peroxidase cycle. Reaction 4: compound I is directly reduced back to the resting state by halides, thereby releasing hypohalous acid. Reactions (1) and (4) constitute the halogenation cycle. Figure adapted from Aldib et al 3 .

Scheme 1 :

 1 Scheme 1: PPIX models serving as a basis for all the other models: a. Minimal model and b. the HRP model.

Figure 3 :

 3 Figure 3: Singly occupied porphyrin orbital (0.02 a. u. isosurface) of Cpd I of the a) MPO and b) MPO imi model, showing a substantial mixing between the he A 1u and A 2u states.

Figure 4 :

 4 Figure 4: a. Quartet and b. doublet transition states of the MPO imi model as representative examples of the TS for Cl -oxidation catalyzed by MPO. The atoms are colored according to the following code: red: oxygen; grey: carbon; white: hydrogen; blue: nitrogen; yellow: sulfur; green: chlorine

Figure 5 :

 5 Figure 5: H-bonding environment involving a conserved Asp and Arg around the proximal His in MPO (C atoms in beige), LPO (pink) and EPO (light blue). The MPO (5FIW) and LPO (2Z5Z) are from crystal structures of respectively human MPO and buffalo LPO. The EPO

Table 1 :

 1 Apparent second order rate constants for Cpd I formation and Cpd I reduction (The MPO, LPO, EPO and HRP data are taken from the more extended Table4of the review of Nicolussi et al.[START_REF] Nicolussi | Posttranslational Modification of Heme in Peroxidases -Impact on Structure and Catalysis[END_REF] , we also point to the references given herein). Values are given at pH = 7.0 and 5.0.

	Substrate	MPO		MPO	MPO	LPO EPO		HRP	CPO 9
				E242Q	M243V
	Cpd I formation (x 10 4 M -1 s -1 )
	H2O2	1800		78	2200	1100 4300		1700	Ca. 400 *
	Cpd I reduction (x 10 4 M -1 s -1 ) by two-electron donors
	pH	7.0 12	5.0 12	7.0	7.0	7.0 13 7.0 14	5.0 14	7.0	4.8
	Chloride	2.5	390	0.0065	-	-	0.31	2.60	-	2.5
	Bromide	110	3000	5.4	13	41	1900	11000	-	38
	Iodide	720	6300	64	1000	1200 9300	> 11000	0.26	n. d.
	Thiocyanate	960	7600	22	10000	12000 10000	> 11000	-	n. d.
	Tyrosine (and nitrite) 77		0.0028	n.d.	11	35		5	n. d.
	-: No activity measured
	n. d.: not determined

Cpd I reduction (x 10 4 M -1 s -1 ) by one-electron donors * Reported in Araiso et al.

Table 4 :

 4 * Estimated from the highest energy scan point obtained in a linear scan of the Cl ➖ --O distance starting from the Cpd I reactant state with bound H2O.Cl -. After identification of the highest energy point, its geometry was subsequently optimized using a frozen optimization procedure to the ground state in which the Cl ➖ --O distance was kept frozen to the one of the highest energy point. Mulliken spin density on Fe, oxygen (O) and on the porphyrin ring (Por) of the oxoiron in the Cl -bound Cpd I reactant complex (Cpd I.H2O.Cl ➖ , quartet spin state).

	Reference state (2) *
	Cpd I + H2O.Cl ➖	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00
	Reaction state ($)
	Cpd I.H2O.Cl ➖	9.35	12.05	19.60	1.42	8.99	2.22	11.02	10.77	21.68
	TS for Cl -oxidation ($)
	Cpd I.H2O.Cl ➖ (TS)	23.36	28.13	37.33	16.37	24.89	13.98	25.42	15.42	44.52
	Product (2)
	Fe(III)+ H2O.ClO ➖	2.03	7.41	n.d.	-2.91	13.47	-5.84	n.d.	n. d.	25.33
	G ‡,$ (with the reactant	14.01	16.08	17.73	14.95	15.90	11.76	14.40	4.65	22.84
	state as reference state)
	b

* The used spin state (see main text and Table

S2

) is given between parenthesis *

Table 5 :

 5 Thermodynamic properties (kcal/mol) calculated for the different models with proximal imidazole, imidazolate and thiolate.It seems to be obvious why MPO is a good oxidant. Its porphyrin ring bears electron withdrawing groups in the form of the ester and sulfonium linkages made to the protein. This resultsin comparison to other imidazole bearing proteins and to CPO, having a proximal thiolatein a high/low Mulliken spin density (Table

	E (kcal/mol)	Minimal	HRP	HRP imi	MPO -	LPO	MPO	MPO imi	CPO	CPO envir
		model			E242Q
	Properties of Cpd 0 (2) *
	Reference state
	Cpd 0 + H + + H2OCl ➖	0.00	0.00	n. d.	0.00	0.00	0.00	n. d.	0.00	0.00
	Stability of Cpd 0
	Cpd I + H2O + H2O.Cl ➖ -288.66	-295.48	n. d.	-283.03	-292.64	-275.81	n. d.	-317.22	-312.27
	Properties of Cpd I (4)
	Reference state
	Cpd I + H2O + H2O.Cl ➖ 0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00
	1-electron affinity Cpd I
	Cpd II + H2O + H2O.Cl -	-123.67	-116.84	-103.19	-126.92	-119.64	-130.83	-116.99	-92.08	-100.33
	Uptake H atom by Cpd I
	Cpd II-OH + H2O +	n. a.	-406.15	-412.85	-405.9	-408.24	-408.30	-414.87	-413.01	-415.65
	H2O.Cl ➖ -H +
	Properties of Cpd I in
	presence of Cl ➖
	Reference state
	Cpd I.H2O.Cl ➖	0.00	0.00	n. d.	0.00	0.00	0.00	n. d.	0.00	0.00
	1 electron affinity Cpd I

* The used spin state (see main text and TableS2) is given between parenthesis n. a.: not available, no convergence during geometry optimization of the Cpd II-OH structure n. a. * : not available, no convergence during geometry optimization of the structure of Cpd II with H2O.Cl ➖ bound n. d.: not determined

Acknowledgement

Computational resources were provided by the Centre de Ressources Informatiques (CRI) of the Université de Lille (France).

Supporting Information

Coordinates of optimized geometries of the redox intermediates: provided as Table S1 in excel-sheet Word document with: Table S2 -Table S5; Text S1 -S2; Figure S1