compare with Hodgson

analyses/Makefile

@@ -8,5 +8,6 @@ all:
 	python plot_DH.py -o ../FigureS1.pdf
 	python plot_CIP.py -o ../FigureS2.pdf
 	python plot_pairs.py -o ../FigureS3.pdf
-	python plot_pot_hammet.py -o ../FigureS4.pdf
-	python plot_cplx_Gx1.py -o ../FigureS5.pdf
+	python plot_pot_vs_Hodgson.py -o ../FigureS4.pdf
+	python plot_pot_hammet.py -o ../FigureS5.pdf
+	python plot_cplx_Gx1.py -o ../FigureS6.pdf

analyses/plot_pot_vs_Hodgson.py

@@ -0,0 +1,109 @@
+import pandas
+import matplotlib.pyplot as plt
+import numpy
+import sys
+import pathlib
+import argparse
+from scipy.spatial import distance_matrix
+from matplotlib.patches import Ellipse
+from nitroxides.commons import dG_DH, AU_TO_ANG, LabelPositioner, AU_TO_EV, EPSILON_R
+from nitroxides.tex import format_longtable
+
+SOLVENT = 'water'
+
+LABELS = {'E_ox': [], 'E_red': []}
+POINTS_POSITION ={'E_ox': [], 'E_red': []}
+LABELS_KWARGS = {'E_ox': [], 'E_red': []}
+LABELS_PATH = {'E_ox': pathlib.Path('pot_hodgson_ox.pos'), 'E_red': pathlib.Path('pot_hodgson_red.pos')}
+
+TO_SHE = 4.4
+
+def plot_family(ax, data: pandas.DataFrame, data_hog: pandas.DataFrame, family: str, column: str, color: str):
+    subdata = data[(data['family'] == family) & (data['solvent'] == SOLVENT)]
+
+    subdata.insert(1, 'compound', [int(n.replace('mol_', '')) for n in subdata['name']])
+    subdata = subdata.join(data_hog[data_hog['E_red'].notnull()].set_index('compound'), on='compound', lsuffix='our', rsuffix='them', how='inner')
+
+    ax.plot(subdata[column + 'our'], subdata[column + 'them'] / 1000 + TO_SHE, 'o', color=color, label=family.replace('Family.', ''))
+
+    for name, e1, e2 in zip(subdata['name'], subdata[column + 'our'], subdata[column + 'them'] / 1000 + TO_SHE):
+        name = name.replace('mol_', '')
+        LABELS[column].append(name)
+        POINTS_POSITION[column].append((e1, e2))
+        LABELS_KWARGS[column].append(dict(color=color, ha='center', va='center'))
+
+def mark_diff(ax, data: pandas.DataFrame, data_hog: pandas.DataFrame, column: str, pos: tuple):
+    subdata = data[data['solvent'] == SOLVENT]
+
+    subdata.insert(1, 'compound', [int(n.replace('mol_', '')) for n in subdata['name']])
+    subdata = subdata.join(data_hog[data_hog['E_red'].notnull()].set_index('compound'), on='compound', lsuffix='our', rsuffix='them', how='inner')
+
+    diff_our_them = subdata[column + 'them'] / 1000 + TO_SHE - subdata[column + 'our']
+
+    ax.text(*pos, '$\\Delta E^0$ = {:.2f} ± {:.2f} V'.format(numpy.mean(diff_our_them), numpy.std(diff_our_them)))
+
+parser = argparse.ArgumentParser()
+parser.add_argument('-i', '--input', default='../data/Data_pot.csv')
+parser.add_argument('-i2', '--input2', default='../data/Data_pot_Hodgson.csv')
+parser.add_argument('-r', '--reposition-labels', action='store_true')
+parser.add_argument('-o', '--output', default='Data_pot_vs_Hodgson.pdf')
+
+args = parser.parse_args()
+
+data = pandas.read_csv(args.input)
+data_hog = pandas.read_csv(args.input2)
+
+figure = plt.figure(figsize=(6, 10))
+ax1, ax2 = figure.subplots(2, 1)
+
+plot_family(ax1, data, data_hog, 'Family.AMO', 'E_ox', 'tab:pink')
+plot_family(ax1, data, data_hog, 'Family.P6O', 'E_ox', 'tab:blue')
+plot_family(ax1, data, data_hog, 'Family.P5O', 'E_ox', 'black')
+plot_family(ax1, data, data_hog, 'Family.IIO', 'E_ox', 'tab:green')
+plot_family(ax1, data, data_hog, 'Family.APO', 'E_ox', 'tab:red')
+
+mark_diff(ax1, data, data_hog, 'E_ox', (5, 5.4))
+
+ax1.set_xlabel('$E^0_{abs}(N^+|N^\\bullet)$ from this work (V)') 
+ax1.set_ylabel('$E^0_{abs}(N^+|N^\\bullet)$ from Hodgson et al. (V)') 
+
+positioner = LabelPositioner.from_file(
+    LABELS_PATH['E_ox'], 
+    numpy.array(POINTS_POSITION['E_ox']), 
+    LABELS['E_ox'], 
+    labels_kwargs=LABELS_KWARGS['E_ox']
+)
+
+if args.reposition_labels:
+    positioner.optimize(dx=1e-3, beta=1e5, krep=1, kspring=1000, c=0.05, b0=0.015)
+    positioner.save(LABELS_PATH['E_ox'])
+
+positioner.add_labels(ax1)
+
+plot_family(ax2, data, data_hog, 'Family.AMO', 'E_red', 'tab:pink')
+plot_family(ax2, data, data_hog, 'Family.P6O', 'E_red', 'tab:blue')
+plot_family(ax2, data, data_hog, 'Family.P5O', 'E_red', 'black')
+plot_family(ax2, data, data_hog, 'Family.IIO', 'E_red', 'tab:green')
+plot_family(ax2, data, data_hog, 'Family.APO', 'E_red', 'tab:red')
+
+mark_diff(ax2, data, data_hog, 'E_red', (2.9, 2.2))
+
+ax2.set_xlabel('$E^0_{abs}(N^\\bullet|N^-)$ from this work (V)') 
+ax2.set_ylabel('$E^0_{abs}(N^\\bullet|N^-)$ from Hodgson et al. (V)') 
+
+positioner = LabelPositioner.from_file(
+    LABELS_PATH['E_red'], 
+    numpy.array(POINTS_POSITION['E_red']), 
+    LABELS['E_red'], 
+    labels_kwargs=LABELS_KWARGS['E_red']
+)
+
+if args.reposition_labels:
+    positioner.optimize(dx=1e-3, beta=1e3, krep=1, kspring=1000, c=0.3, b0=0.05, scale=[3, 1])
+    positioner.save(LABELS_PATH['E_red'])
+
+positioner.add_labels(ax2)
+
+plt.tight_layout()
+plt.legend()
+figure.savefig(args.output)

analyses/pot_hodgson_ox.pos

@@ -0,0 +1,54 @@
+,label,x,y
+0,1,4.984310763602672,5.103258287011272
+1,2,5.096418431400864,5.18509275270964
+2,3,5.17227559718236,5.209345074989223
+3,4,5.13133209490243,5.088688196291432
+4,5,5.164241262338934,5.17668475202551
+5,6,5.1965152780896116,5.1464163245031935
+6,7,5.081765388515453,5.129742064940129
+7,8,5.186053059826803,5.265569956730343
+8,9,5.2825533183131235,5.322414317797528
+9,10,5.312828850101342,5.4303624285555365
+10,11,5.304195872991947,5.260848118327276
+11,12,5.314434422329497,5.343011521261141
+12,13,5.135469387414481,5.20786220920241
+13,14,5.068756307035109,5.163685729928762
+14,15,5.253050082343705,5.295081507629994
+15,16,5.15217674163271,5.146095246970255
+16,17,5.11358940617911,5.118295783664926
+17,18,5.164664659011584,5.110749890088508
+18,19,5.128504099077829,5.171937179359233
+19,20,5.330956338280313,5.382044042798034
+20,22,5.337155121777711,5.412001552893007
+21,23,5.255069739374381,5.402274320959204
+22,24,5.278854139406711,5.461402116515664
+23,25,5.226120738896293,5.327140906041079
+24,26,5.214981951821337,5.36421145785555
+25,27,5.235854037303373,5.382414308430322
+26,28,5.228015117079243,5.4350171051657785
+27,29,5.275790551552028,5.428560304208335
+28,30,5.376606486630395,5.413484971394564
+29,31,5.3095102886854075,5.4620871656647685
+30,32,5.34505825931232,5.526042609490015
+31,33,5.3562292352369605,5.450194255911371
+32,34,5.285638590307379,5.371601931757725
+33,35,5.341433322070314,5.48549285160866
+34,36,5.229290622554556,4.869238268074303
+35,37,5.248875903305158,4.892030116981909
+36,38,5.313744718097044,4.922994362194325
+37,39,5.308361299192163,4.904766528182234
+38,40,5.299100548466084,4.872824892856739
+39,41,5.2798126085175605,4.938697661003501
+40,42,5.331222409705173,4.888784921023083
+41,43,5.300572538261673,4.966950636776302
+42,44,5.274362690958083,4.907617807315995
+43,45,5.310935529193714,5.033601511436155
+44,46,5.328167294681219,5.002083765806227
+45,47,5.273660195232492,5.004679134000951
+46,48,5.324060356015427,4.949508745798863
+47,49,5.134313541156714,4.780915620663185
+48,50,5.236362835855231,4.824887399091428
+49,51,5.1727792047245345,4.757134557290438
+50,52,5.361624613447284,4.989073325652078
+51,53,5.377408668428593,4.948788306227779
+52,54,5.346681130885673,4.918870874722461

analyses/pot_hodgson_red.pos

@@ -0,0 +1,54 @@
+,label,x,y
+0,1,2.861830371021398,2.87238785691857
+1,2,2.9462867550838117,3.02202051767794
+2,3,3.012593200919952,3.1163648302605145
+3,4,2.9735585761839607,3.0096828105213245
+4,5,2.994645238661583,3.0549150927534594
+5,6,3.0048023065443736,2.981817680242245
+6,7,2.788218155772184,2.6012347315320508
+7,8,3.0400312391692563,3.1132311675379674
+8,9,3.113522249355996,3.187575696506334
+9,10,3.1661616704100704,3.1966179458007016
+10,11,3.119639055038228,3.091968425403241
+11,12,3.139031559840304,3.222477838122566
+12,13,2.98117663280009,2.9343447905143845
+13,14,2.835614423250416,2.728587871201881
+14,15,2.927310127829675,2.9229691414756416
+15,16,2.8822477439175036,2.7604625156265006
+16,17,2.8963195706568206,2.866823298275264
+17,18,2.9358696808758102,2.7549862157719605
+18,19,2.739988596171574,2.6839248281411865
+19,20,3.023606549001566,3.045187676744988
+20,22,3.1053762252230306,2.9715233699701855
+21,23,2.9537904897271186,2.9494950501409622
+22,24,3.0878473959078563,2.8621086226537904
+23,25,2.941695339032419,2.8467337698801543
+24,26,2.965368898260228,2.7948029282625084
+25,27,2.9966255150977665,2.789624306533792
+26,28,2.7776454307857867,2.45819624597765
+27,29,3.0650337212172736,2.7456282886375263
+28,30,3.1062964656053107,2.0324516522139806
+29,31,3.1087822567244237,2.665181382035088
+30,32,3.0957070067112467,2.751124186583203
+31,33,3.038908087085751,2.920016434227756
+32,34,3.0156221249302617,2.859115663767425
+33,35,3.0465752432361066,2.6480985546123312
+34,36,3.0938711600146793,3.13159476050939
+35,37,3.1429444541601166,3.1394699680604123
+36,38,3.075326920993851,2.9543899480580014
+37,39,3.164464694009839,3.001286303781044
+38,40,3.2001951712924663,3.1023961632764516
+39,41,3.228560540127451,2.9244207879305306
+40,42,3.0806733494477903,2.559037183165133
+41,43,3.1386570153738784,2.823334289168122
+42,44,3.2173350409594157,3.22256074968582
+43,45,3.208997871425079,3.016403312953077
+44,46,3.1634105470333127,2.8849620869816475
+45,47,3.135948910259371,2.9661581900124285
+46,48,3.17087091123991,3.099059659918973
+47,49,3.066145860765712,3.1413341189937607
+48,50,3.094798980004314,3.0468976331547175
+49,51,3.0491641724956624,3.0010284612901614
+50,52,3.1964324553744383,2.903929017675193
+51,53,3.145117683878293,2.1676902317375712
+52,54,3.1571536897719925,2.0184152910141737

data/Data_pot_Hodgson.csv

@@ -0,0 +1,55 @@
+compound,E_ox,E_red
+1,696.5,-1552.9
+2,762.1,-1482.5
+3,793.7,-1387.6
+4,703.7,-1455.3
+5,779.8,-1436.3
+6,756.4,-1413.4
+7,741.3,-1783.3
+8,878.5,-1347.9
+9,914,-1270
+10,1016.1,-1217.1
+11,853.8,-1315.2
+12,946.1,-1255.1
+13,793.9,-1479.52
+14,763.7,-1621.3
+15,909.9,-1504.8
+16,729.5,-1591.4
+17,726.9,-1552.6
+18,718.5,-1598
+19,765.3,-1722.7
+20,997.2,-1396
+21,1043.9,
+22,1023,-1349.6
+23,1000.3,-1504.8
+24,1055.7,-1536.8
+25,941.7,-1525.7
+26,974.4,-1554
+27,971,-1561
+28,1024.2,-1935.5
+29,1016.2,-1667.2
+30,1022.9,-2331
+31,1072.2,-1754.7
+32,1106.5,-1616.7
+33,1056.6,-1423.7
+34,965.2,-1497.5
+35,1100.3,-1757.6
+36,469.1,-1336
+37,481,-1313.4
+38,529.8,-1390.5
+39,512.7,-1355.6
+40,486.1,-1328.4
+41,527.2,-1475.4
+42,510.8,-1792
+43,554.2,-1542.6
+44,509.6,-1226.3
+45,619.2,-1402
+46,593.4,-1468.2
+47,601.1,-1384.2
+48,540.7,-1328.7
+49,368.7,-1317.9
+50,439.9,-1346.7
+51,359.2,-1369
+52,575.2,-1449.6
+53,555,-2236.5
+54,514,-2343.1

nitroxides.tex

@@ -376,7 +376,7 @@ \section{Methodology} \label{sec:methodo}
 \label{fig:nitroxides}
 \end{figure}
 
-Geometry optimizations and subsequent vibrational frequency calculations were performed at the $\omega$B97X-D/6-311+G(d) level in water and acetonitrile (described using the SMD \cite{marenichUniversalSolvationModel2009} approach) with Gaussian 16 C02 \cite{g16}. With other possible candidates, this functional have been demonstrated to provide reliable results  \cite{flores-leonarFurtherInsightsDFT2017,maierG4AccuracyDFT2020}.  For compound \textbf{1}-\textbf{54}, the geometries obtained by Hodgson et al. \cite{hodgsonOneElectronOxidationReduction2007} have been used as a starting point, taking advantage of their extensive conformational search. All radical forms are considered to have a doublet ground state  [$\braket{S^2}=\frac{3}{4}$]. Then, the same calculations were preformed in acetonitrile for the subset of compounds for which experimental redox potentials are available (listed in Fig.~\ref{fig:nitroxides}). Finally, to study the influence of the substituent on the redox potential with the model presented in Section \ref{sec:eleczhang}, single point calculation are performed at the $\omega$B97X-D/6-311+G(d) level in gas phase, using the optimized geometries of  the radical states of each nitroxides (in water) in which $>$\ce{N-O^.} moiety is substituted by \ce{CH_2} (the rest of the geometry is kept fixed).
+Geometry optimizations and subsequent vibrational frequency calculations were performed at the $\omega$B97X-D/6-311+G(d) level in water and acetonitrile (described using the SMD \cite{marenichUniversalSolvationModel2009} approach) with Gaussian 16 C02 \cite{g16}. With other possible candidates, this functional have been demonstrated to provide reliable results  \cite{flores-leonarFurtherInsightsDFT2017,maierG4AccuracyDFT2020} (see also Fig.~S4).  For compound \textbf{1}-\textbf{54}, the geometries obtained by Hodgson et al. \cite{hodgsonOneElectronOxidationReduction2007} have been used as a starting point, taking advantage of their extensive conformational search. All radical forms are considered to have a doublet ground state  [$\braket{S^2}=\frac{3}{4}$]. Then, the same calculations were preformed in acetonitrile for the subset of compounds for which experimental redox potentials are available (listed in Fig.~\ref{fig:nitroxides}). Finally, to study the influence of the substituent on the redox potential with the model presented in Section \ref{sec:eleczhang}, single point calculation are performed at the $\omega$B97X-D/6-311+G(d) level in gas phase, using the optimized geometries of  the radical states of each nitroxides (in water) in which $>$\ce{N-O^.} moiety is substituted by \ce{CH_2} (the rest of the geometry is kept fixed).
 
 Since all thermochemical quantities are $\kappa$-dependent, analyses were performed using custom Python scripts. When required (e.g., in Eq.~\eqref{eq:dh}), the value of $a$ (the radius of the solute cavity) is taken as half the largest distance between two atoms in the molecule. Although this is an approximation, it provides a consistent method to treat all molecules proportionally to their size and is consistent with other publications \cite{matsuiDensityFunctionalTheory2013}.  Furthermore, a value of $\varepsilon_{r,wa}=80$ for water and $\varepsilon_{r,ac}=35$ for acetonitrile is used. These relative permittivities correspond to those of the pure solvents and are known to be lower in the respective electrolyte solutions \cite{silvaTrueHuckelEquation2022}. These variations can be substantial; for example, $\varepsilon_r \approx 70$ for a solution containing \SI{1}{\mol\per\kilo\gram} of \ce{NaCl} in water \cite{kontogeorgisDebyeHuckelTheoryIts2018, silvaTrueHuckelEquation2022}, but they depend on the nature of the electrolyte, so it was not considered here.
 
@@ -410,7 +410,7 @@ \subsection{Structure-activity relationships} \label{sec:sar}
 As a consequence, \textbf{55} exhibits the highest oxidation and reduction potentials among all the compounds studied in this paper.
 
 
-To elucidate these effects, attempts were made to correlate both potentials with Hammett constants for P5O and P6O, but the correlations were found to be very weak, especially for reduction (see Fig.~S4). The electrostatic interaction model [Eq.~\eqref{eq:Er}] provides more insights. Results are presented in Fig.~\ref{fig:corr} (see also Table S3). However, this model fails to account for the effect of substituting methyl groups with ethyl groups. Moreover, including the disubstituted compounds (e.g., \textbf{9}) worsens the correlation ($R^2 \sim 0.5$ and 0.3 for oxidation and reduction, respectively). Compounds \textbf{56} and \textbf{58} remain outliers for reduction. Therefore, all three sets of compounds were treated as outliers.
+To elucidate these effects, attempts were made to correlate both potentials with Hammett constants for P5O and P6O, but the correlations were found to be very weak, especially for reduction (see Fig.~S5). The electrostatic interaction model [Eq.~\eqref{eq:Er}] provides more insights. Results are presented in Fig.~\ref{fig:corr} (see also Table S3). However, this model fails to account for the effect of substituting methyl groups with ethyl groups. Moreover, including the disubstituted compounds (e.g., \textbf{9}) worsens the correlation ($R^2 \sim 0.5$ and 0.3 for oxidation and reduction, respectively). Compounds \textbf{56} and \textbf{58} remain outliers for reduction. Therefore, all three sets of compounds were treated as outliers.
 On the positive side, this model helps explain some of the effects mentioned above: the increase in oxidation (and reduction) potential for aromatic compounds correlates with an increase in quadrupole moment ($Q_{xx} > \SI{5}{\elementarycharge\bohr\squared}$ for most member of IIO or APO), while the modification due to donor/acceptor substituents is linked to changes in the dipole moment. For example, aromatic compounds that have \ce{NH2} has substituent (\textit{e.g.}, \textbf{51}) are characterized by  $\mu_{x} < 0$, which gets larger for coumpounds have \ce{COOH} (\textit{e.g.}, \textbf{39}) or \ce{NO2} (\textit{e.g.}, \textbf{54}). It also accounts for some effects due to the position of the substituent (see, e.g., \textbf{49}-\textbf{51}), which was not the case with the original model by Zhang and co-workers (resulting in weak correlations, $R^2 \leq 0.3$).
 Finally, although it is not directly applicable to charged substituents (\textbf{11}, \textbf{21}, and \textbf{35}), for which the multipole moments are ill-defined, the leading term $q/r$ would result in a positive contribution to $E_r$, which correlates well with the increase in oxidation and reduction potential for these compounds.
 

nitroxides_SI.tex

@@ -37,6 +37,13 @@
 \caption{Evolution of the pairing energy, $\Delta G^\star_{pair}$ (computed with Eq.~(11) of the main text) between two ions, with the ratio between the radii of the two ions ($\chi$) and for increasing dipole cavity shapes ($s_2$). Their charge is set to 1 and -1, and two possible scaling for  the close contact distance is used ($s_1$, left and right). The equilibrium constant, $K_{pair} = e^{-\frac{\Delta G^\star_{pair}}{RT}}$, is reported. }
 \end{figure}
 
+
+\begin{figure}[!h]
+	\centering
+	\includegraphics [width=.7\linewidth]{FigureS4}
+	\caption{Comparison between the absolute oxidation (ltop) and reduction (bottom) potentials obtained by Hodgson \emph{et al.} (at the G3(MP2)-RAD:MP2/6-311+G(3df,2p) level, using solvation energy computed at the B3LYP/6-31G(d) level in water [CPCM]) \cite{hodgsonOneElectronOxidationReduction2007}, and in this work (at the $\omega$B97X-D/6-311+G(d) level in water [SMD]). The DH correction has \textbf{not} been applied. The large differences (e.g., for the oxidation of APOs) are due to change in geometries after re-optimization.}
+\end{figure}
+
 \clearpage
 \begin{longtblr}[caption={Radii ($a$, in \si{\angstrom}) for all oxidized states of compounds \textbf{1}-\textbf{61} and corresponding absolute redox potentials ($E^0_{abs}$, in \si{\volt}), as computed at the $\omega$B97X-D/6-311+G(d) level in water (SMD), with $[\ce{X}]=\SI{0}{\mole\per\liter}$.}]{colspec={>{\bfseries}lX[c]X[c]X[c]X[c]X[c]},width = 0.85\linewidth, rowhead=1}
 	\hline
@@ -150,7 +157,7 @@
 		\ce{NO2} & 0.78 & 0.71 \\
 		\hline
 	\end{tblr}
-	\includegraphics [width=\linewidth]{FigureS4}
+	\includegraphics [width=\linewidth]{FigureS5}
 	\caption{Correlation between absolute oxidation (left) and reduction (right) and the Hammet constant of their substituent (para, $\sigma_p$, or meta, $\sigma_m$, from Ref.~\cite{hanschSurveyHammettSubstituent1991} and reported in the table) for compounds of the P5O (black markers, correlated with $\sigma_m$) and P6O (blue markers, with correlated with $\sigma_p$) families.}
 \end{figure}
 
@@ -317,7 +324,7 @@
 
 \begin{figure}[!h]
 \centering
-\includegraphics [width=\linewidth]{FigureS5}
+\includegraphics [width=\linewidth]{FigureS6}
 \caption{Value of the complexation free Gibs energy change for $K_{01}$ (round markers, $\bullet$) and $K_{21}$ (square markers, $\blacksquare$), as computed at the $\omega$B97X-D/6-311+G(d) level in water (top) and acetonitrile (bottom) using SMD at two concentration: $[X]=\SI{1}{\mole\per\liter}$ (filled markers) and  $[X]=\SI{0.1}{\mole\per\liter}$ (empty markers). }
 \end{figure}