text:
\section*{Introduction}

Polymer electrolytes are promising candidates for next generation lithium-ion battery technology due to their low cost, safety, and manufacturing compatibility. The major challenge with the current polymer electrolytes is their low ionic conductivity, which limits the usage in real world applications. \cite{hallinan2013polymer,agrawal2008solid,ngai2016review} This limitation has motivated tremendous research efforts to explore new classes of polymers via both experiments \cite{pesko2016effect,tominaga2010alternating,meabe2017polycondensation,hatakeyama2020ai} and atomic scale simulations \cite{webb2015systematic,savoie2017enhancing,france2019effect}. However, the simulation of ionic conductivity is extremely expensive due to the amorphous nature\footnote{Strictly speaking, polymers have different degrees of crystallinity and are not completely amorphous. In this work, we use the term amorphous materials to represent a range of solids from fully amorphous to semi-crystalline. They cannot be represented by a simple periodic structure and thus are generally difficult to simulate.} of polymer electrolytes and the diversity of timescales involved in their dynamics, drastically limiting the ability to employ high throughput computational screening approaches. For instance, recent studies \cite{webb2015systematic,savoie2017enhancing,france2019effect} exploring polymer electrolytes with classical molecular dynamics (MD) only simulated around ten polymers. In contrast, a study that applies machine learning methods to literature data is able to explore a larger chemical space \cite{hatakeyama2020ai}, but it is limited by the diversity of polymers that have been studied in the past. The exploration beyond known chemical spaces would require a significant acceleration of the computational screening of polymer electrolytes.

There are two major reasons for the large computational cost for simulating the ionic conductivity of polymer electrolytes with MD. First, the amorphous structure of polymer electrolytes can only be sampled from a random distribution using, e.g., Monte Carlo algorithms, and yet this initial structure has a significant impact on the simulated ionic conductivity due to the lack of ergodicity in the MD simulation \cite{molinari2018effect,france2019effect}. Multiple simulations starting from independent configurations are therefore required in order to properly sample the phase space and reduce statistical noise. Second, the slow relaxation of polymers requires long MD simulation time to achieve convergence for ionic conductivity (on the orders of 10's to 100's of ns), so each MD simulation is also computationally expensive \cite{webb2015systematic,france2019effect}.

Machine learning (ML) techniques have been widely used to accelerate the screening of ordered materials \cite{butler2018machine, schmidt2019recent}, but most previous studies implicitly \cite{gomez2016design,ye2018deep,ahmad2018machine,de2016statistical} assume that the properties used to train the ML models are generated through a deterministic, unbiased process. However, the MD simulation of complex materials like amorphous polymers is intrinsically stochastic, and obtaining data with low statistical uncertainties by running repetitive simulations is impractical at a large scale due to the large computational cost. An alternative approach is to reduce the accuracy requirements for individual MD simulations and learn to reduce the random and systematic errors with large quantities of less expensive, yet imperfect data. It has previously been demonstrated that ML models can learn from noisy data and recover the true labels for images \cite{rolnick2017deep} and graphs \cite{du2018robust}. Past works have also shown that systematic differences between datasets can be learned employing transfer learning techniques. \cite{yamada2019predicting,jha2019enhancing,smith2019approaching,wu2019machine} Inspired by these results, we hope to significantly reduce the computational cost for simulating the transport behavior of polymers by adopting a noisy, biased simulation scheme with short, unconverged MD simulations.

In this work, we aim to accelerate the high throughput computational screening of polymer electrolytes by learning from a large amount of biased, noisy data and a small number of unbiased data from molecular dynamics simulations. Despite the large random errors caused by the dependence on initial structure, we only perform one MD simulation for each polymer, and learn a shared model across polymers to reduce the random error and recover true properties that one would obtain from repetitive simulations. To reduce the long MD simulation time, we perform large quantities of short, unconverged MD simulations and a small number of long, converged simulations. We then employ multi-task learning to learn a correction from the short simulation properties to long simulation properties. We find that our model achieves a prediction error with respect to the true properties smaller than the random error from a single MD simulation, and it also corrects the systematic errors from unconverged simulations better than a linear correction. Combining the reduction of both random and systematic errors, we successfully screen a space of 6247 polymers and discover the best polymer electrolytes from the space, which corresponds to a 22.8-fold acceleration compared with simulating each polymer directly with one long simulation. Finally, we extrapolate our model into a larger space of 53362 polymers with more complex monomer structures, and we validate our predictions with 31 experimentally studied polymers.

\section*{Results}

\begin{figure*}[tb]
  \centering
  \includegraphics[width=\linewidth]{illustrative.pdf}
  \caption{Illustration of the polymer space and the learning framework. (a) The polymer space and molecular dynamics simulation workflow. (b) Ionic conductivity as a function as simulation time from 6 independent 5 ns MD runs for the same polymer, showing the random errors caused by the amorphous structure initialization. (c) Ionic conductivity as a function as simulation time for 5 different polymers, showing the long simulation time needed for convergence. (Polymer structures for (b, c) are provided in the supplementary information.) (d) Multi-task learning framework to reduce the random and systematic errors from MD simulations.}
  \label{fig:illustrative}
\end{figure*}

\textbf{Polymer space and sources of errors.} The polymer space we aim to explore is defined in Fig. \ref{fig:illustrative}(a), which considers both the synthesizability of polymers and their potential as electrolytes. In general, it is difficult to determine the synthesizability, especially the polymerizability, of an unknown polymer. Here we focus on a well established condensation polymerization route using carbonyl dichloride and comonomers containing any combination of two primary hydroxyl, amino, or thiol groups to form poly-carbonates, ureas, dithiocarbonates, urethanes, thiourethanes, and thiocarbonates.  This scheme does not guarantee polymerizability, but provides a likely route for lab synthesis. The carbonyl structure ensures a minimum capability to solvate Li-ions as an electrolyte, and it also allows for the maximum diversity of polymer backbones. The monomers are sampled from a large pharmaceutical database \cite{irwin2005zinc} to ensure its structures are realistic. After obtaining the molecular structure of the polymer, we sample its 3D amorphous structure with a Monte Carlo algorithm, insert \SI{1.5}{mol} lithium bis(trifluoromethanesulfonyl)imide (LiTFSI) salt per kilogram of polymer, perform a \SI{5}{ns} MD equilibration, and finally run the MD simulation to compute its transport properties like conductivity.

There are mainly two types of errors in this workflow. In the scope of this work, we call random errors the ones that can be eliminated by running repetitive simulations on the same polymer, and systematic errors those that cannot be eliminated. The major source of random error is the sampling of initial amorphous structure of the polymer. In Fig. \ref{fig:illustrative}(b), we show the conductivities computed from 6 different random initializations for the same polymer, which has a large standard deviation of \SI{0.094}{log_{10}(S/cm)} in the log scale at \SI{5}{ns}. This error comes from the lack of ergodicity of MD simulations for polymers -- the large scale amorphous structure of the polymers usually does not change significantly at the timescale that can be achieved with MD. The systematic errors mainly come from the long MD simulation time needed to obtain the converged conductivity. Fig. \ref{fig:illustrative}(c) shows the value of conductivity as a function of the simulation time for 5 different polymers, which slowly converges as the simulation progresses. This slow convergence introduces a systematic error of ionic conductivity with any specified simulation time with respect to the converged conductivity. On average, there is a \SI{0.435}{log_{10}(S/cm)} difference in the log scale between a \SI{5}{ns} and a \SI{50}{ns} simulation for these 5 polymers. Here, we use the \SI{50}{ns} simulation results as the converged values, although it is not fully converged for some polymers. Based on our comparison with respect to experimental values reported in literature \cite{pesko2016universal,zheng2018optimizing,tominaga2017ion,meabe2017polycondensation,mindemark2016allyl,fonseca2006development,itoh2018thermal,pehlivan2010ionic,he2017carbonate,pesko2016effect,doeff2000transport,silva2006novel} in Fig. \ref{fig:SI-md-50ns-exp}, the \SI{50}{ns} simulation has a reasonable agreement except for polymers with very low conductivity. Note that even 50 ns conductivities have large random errors similar to the 5 ns conductivities, since the random errors are mainly caused by the large scale amorphous structures that do not change significantly with long simulation time.

 \textbf{Multi-task model to reduce errors.} These two types of errors introduce significant computational costs to achieve an accurate calculation of ionic conductivity for individual polymers, because such a calculation requires repetitive simulations on the same polymer that are also individually expensive. Here we attempt to reduce these errors by learning a shared model across the polymer space. To achieve this goal, we develop a multi-task graph neural network architecture (Fig. \ref{fig:illustrative}(d)) to learn to reduce both random and systematic errors from MD simulations. We first encode the monomer structure as a graph $\mathcal{G}$ (details of the encoding discussed in the methods section) and use a graph neural network $G$ to learn a representation for the corresponding polymer, $
 	\bm{v}_{\mathcal{G}} = G(\mathcal{G})$. Here we use a CGCNN \cite{xie2018crystal} as $G$, similar to previous works that employ graph convolutional networks (GCNs) in polymers \cite{zeng2018graph,st2019message}.

 To build a predictor that reduces random errors, we use the robustness of neural networks against random noises in the training data, previously demonstrated in images \cite{rolnick2017deep} and graphs \cite{du2018robust}. We assume that there exists a true target property (e.g. conductivity) that is uniquely determined by the structure of the polymer (which would require infinite repetitive simulations to obtain), and the computed target property from MD is slightly different from the true property due to the random errors in the simulation. This assumption can be written as,
 \begin{equation}\label{eq:noise}
 	t = f(\mathcal{G}) + \epsilon,
 \end{equation}
where $t$ is the target property computed from MD, $f$ is a deterministic function mapping from monomer structure to true polymer property, and $\epsilon$ is a random variable independent of $\mathcal{G}$ with zero mean. Note that $\epsilon$ should be a function of $\mathcal{G}$ in principle, but similar noise is observed across polymers as shown in Fig. \ref{fig:SI-true-error} and assuming $\epsilon$ is independent of $\mathcal{G}$ simplifies our analysis. By regressing over $t$, it is possible to learn $f(\mathcal{G})$ even when the noise is large \cite{rolnick2017deep} if enough  training data is available. To generate a large amount of training data, since \SI{50}{ns} simulations are too expensive practically, we use less accurate \SI{5}{ns} simulations to generate training data and use a network $g_1$ to predict $t_{\SI{5}{ns}}$ with the graph representation,
\begin{equation}\label{eq:network-5ns}
	y_{\SI{5}{ns}} = g_1(\bm{v}_{\mathcal{G}}).
\end{equation}

With enough training data generated using the affordable $\SI{5}{ns}$ simulations, we can learn an approximation to the true property function $f_{\SI{5}{ns}}$ despite the random errors. However, there is a systematic error between $f_{\SI{5}{ns}}$ and $f_{\SI{50}{ns}}$ due to the slow relaxation of polymers. To correct this error, we perform a small amount of \SI{50}{ns} simulation to generate data for the converged conductivities. This correction can then be learned with a linear layer $g_2$ using both predictions from $\SI{5}{ns}$ simulations and the graph representations,
\begin{equation}\label{eq:network-50ns}
	y_{\SI{50}{ns}} = g_2(\bm{v}_{\mathcal{G}} \mathbin\Vert y_{\SI{5}{ns}}),
\end{equation}
where $\mathbin\Vert$ denotes concatenation.

Finally, the two datasets, a larger $\SI{5}{ns}$ dataset and a smaller $\SI{50}{ns}$ dataset, can be trained jointly using a combined loss function,
\begin{equation}
	\mathrm{Loss} = (1 - w) \cdot \frac{1}{N_{\SI{5}{ns}}}\sum_{\mathcal{G}_{\SI{5}{ns}}}(y_{\SI{5}{ns}} - t_{\SI{5}{ns}})^2 +
		w \cdot \frac{1}{N_{\SI{50}{ns}}}\sum_{\mathcal{G}_{\SI{50}{ns}}}(y_{\SI{50}{ns}} - t_{\SI{50}{ns}})^2,
\end{equation}
where $w$ is a weight between 0 and 1.

Using an iterative scheme, we sampled the entire polymer space in Fig. \ref{fig:illustrative}(a) with both $\SI{5}{ns}$ and $\SI{50}{ns}$ simulations. The $\SI{5}{ns}$ dataset includes 876 polymers and the $\SI{50}{ns}$ dataset includes 117 polymers. Note that we only simulate each polymer once so there is no duplicate in both datasets. We leave 10\% of the polymers in both datasets as test data, and use 10-fold cross validation on the rest of the data to train our models. Due to the small size of the $\SI{50}{ns}$ dataset, we use stratified split while dividing the data to ensure that the training, validation, test data contain polymers with the full range of conductivities \cite{wu2018moleculenet}. In the next sections, we first demonstrate the performance of our model based on these two datasets and then discuss the iterative screening of the polymer space.

\begin{figure*}[tb]
  \centering
  \includegraphics[width=\linewidth]{random_error.pdf}
  \caption{Performance on reducing random errors. (a) Mean absolute errors (MAEs) on a toy dataset to predict LogP with increasing noises in training data. Blue line denotes MAEs with respect to true LogP values, green line denotes MAEs with respect to noisy LogP values, and dashed line denotes the mean absolute deviation (MAD) of the Gaussian noise. (b) Scatter plot comparing the predicted conductivity and computed conductivity from one initialization (config A) in the test dataset. (c) Scatter plot comparing conductivities from two independent initializations (config A and config B) for the same polymers in the test dataset.}
  \label{fig:random-error}
\end{figure*}

\textbf{Performance on reducing random errors.} To demonstrate that our model can recover the true properties from noisy data, we first study a toy dataset for which we have access to the true property $f(\mathcal{G})$ in Eq. \ref{eq:noise}. We use the same dataset from $\SI{5}{ns}$ simulations and compute the partition coefficient, LogP, of each polymer using Crippen's approach \cite{wildman1999prediction, rdkit}, which uses an empirical equation whose output is fully determined by the molecular structure. Then, we add different levels of Gaussian random noise into the LogP values to imitate the random errors in simulated conductivities. Here, we only use the $g_1$ branch of our model, i.e. $w=0$, to predict LogP values from the synthesized noisy data. Fig.~\ref{fig:random-error}(a) shows the true mean absolute errors (MAEs) with respect to the original LogP values and apparent MAEs with respect to the noisy LogP values as a function of the standard deviation of the Gaussian noise, on a test dataset including 86 polymers. We observe that the true MAEs become smaller than the mean absolute deviation (MAD) of the Gaussian noise when the noise standard deviation is larger than 0.08. This result shows that our model predicts LogP more accurately than performing a noisy simulation of LogP due to the existence of large random error in the simulation.

We cannot use the same approach to evaluate the model performance on predicting simulated $\SI{5}{ns}$ conductivities because we do not have access to the true conductivities. Instead, we can only compute the apparent MAE on 86 test data, which is \SI{0.117}{log_{10}(S/cm)} in Fig. \ref{fig:random-error}(b). To compare our prediction performance with performing a $\SI{5}{ns}$ MD simulation, we run another independent MD simulation for each test polymer and compare the MAE between the two simulations. The MAE in Fig. \ref{fig:random-error}(c) is \SI{0.132}{log_{10}(S/cm)}, slightly larger than the apparent prediction MAE. This result indicates that our prediction of the noisy conductivity also outperforms an independent MD simulation due to its large random noise, similar to the LogP prediction. To estimate the true prediction performance with respect to the inaccessible true conductivity, we need to assume that the random errors for $\SI{5}{ns}$ MD conductivity follow a Gaussian distribution, which is approximately correct (Fig. \ref{fig:SI-true-error}). We could then estimate the true root mean squared error (RMSE) to be \SI{0.085}{log_{10}(S/cm)}, smaller than the standard deviation of the Gaussian noise \SI{0.117}{log_{10}(S/cm)} (detailed calculations can be found in supplementary note 1). This finding shows that the error in our prediction of \SI{5}{ns} simulated conductivity is smaller than a \SI{5}{ns} MD simulation due to its random errors.

\begin{figure*}[tb]
  \centering
  \includegraphics[width=\linewidth]{systematic_error.pdf}
  \caption{Performance on correcting systematic errors. (a-b) Scatter plots showing the interpolation (a) and extrapolation (b) performance of the model on test data. Blue and green dots present the results of $\SI{5}{ns}$ MD simulations and ML predictions compared with $\SI{50}{ns}$ MD conductivities, respectively. The error bars represent the standard deviations of predictions from 10-fold cross validation. (c) Change of interpolation performance with different number of CV data. The red, blue, green lines denote the MAEs of linear correction, multi-task model, and single-task model predicting $\SI{50}{ns}$ conductivity. The error bars represent the standard deviations of MAEs from 10-fold cross validation.}
  \label{fig:systematic-error}
\end{figure*}

\textbf{Performance on correcting systematic errors.} In addition to reducing random errors, our model is also able to learn the systematic difference between $\SI{5}{ns}$ and $\SI{50}{ns}$ MD simulated properties with the multi-task scheme. After co-training our model with both $\SI{5}{ns}$ and $\SI{50}{ns}$ datasets, we present the predictions on 11 test data from $\SI{50}{ns}$ MD in Fig. \ref{fig:systematic-error}(a). Compared with the original $\SI{5}{ns}$ conductivities, our model corrects the systematic error and achieves a MAE of \SI{0.076}{log_{10}(S/cm)} by averaging the predictions from 10-fold cross validations. It is clear that the model corrects the systematic error by learning a customized correction to each polymer, which is better than an overall linear correction which gives a MAE of \SI{0.152}{log_{10}(S/cm)}. Note that this MAE does not include random errors, because our $\SI{5}{ns}$ and $\SI{50}{ns}$ conductivities are computed from the same random initial structures. The results in Fig. \ref{fig:systematic-error}(a) represent the interpolation performance of our model since we randomly split our data. To further study the extrapolation performance, we perform the same co-training but reserve the top 10 polymers with highest conductivity as test data. In Fig. \ref{fig:systematic-error}(b), we find that by training with low-conductivity polymers, the model underestimates the $\SI{50}{ns}$ conductivity and achieves a MAE of \SI{0.182}{log_{10}(S/cm)}. This underestimation is due to the larger systematic error between $\SI{50}{ns}$ and $\SI{5}{ns}$ conductivities in training data, caused by slow relaxations in low-conductivity polymers. Nevertheless, the model still performs better than a linear correction that only has access to the training data, which has a MAE of \SI{0.275}{log_{10}(S/cm)}.

\begin{table*
\caption{Comparison of the mean absolute errors (MAEs) on predicting $\SI{50}{ns}$ MD simulated properties between different approaches. For each property, interpolation and extrapolation performance are represented by labels without and with the ${}^*$ symbol. Uncertainties are the standard deviations of MAEs from 10-fold cross validation (CV).\label{tab:comparison}}
\begin{ruledtabular}
\begin{tabular}{p{2.5cm}*{8}{p{1.2cm}}}

  Method & $\sigma$ & $\sigma^*$ & $D_{\mathrm{Li}}$ & $D_{\mathrm{Li}}^*$ & $D_{\mathrm{TFSI}}$ & $D_{\mathrm{TFSI}}^*$ & $D_{\mathrm{Poly}}$ & $D_{\mathrm{Poly}}^*$ \\
  5 ns & 0.528 & 0.278 & 0.503 & 0.419 & 0.455 & 0.249 & 0.612 & 0.528\\
  5 ns (linear) & 0.152 & 0.275 & 0.148 & 0.247 & 0.096 & 0.297 & 0.072 & 0.110\\
  ML CV & $0.093 \pm 0.017$ & $0.186 \pm 0.053$ & $0.106 \pm 0.016$ & $0.209 \pm 0.050$ & $0.101 \pm 0.020$ & $0.181 \pm 0.028$ & $0.072 \pm 0.019$ & $0.114 \pm 0.030$\\
  ML average & 0.076 & 0.182 & 0.080 & 0.202 & 0.075 & 0.171 & 0.056 & 0.104\\

\end{tabular}
\end{ruledtabular}
\end{table*}

In Table \ref{tab:comparison}, we study how the systematic error correction performs for other transport properties, including lithium ion diffusivity ($D_{\mathrm{Li}}$), TFSI diffusivity ($D_{\mathrm{TFSI}}$), and polymer diffusivity ($D_{\mathrm{Poly}}$). Both interpolation and extrapolation performances are reported similar to the results of conductivity. To better evaluate the uncertainties caused by the small $\SI{50}{ns}$ dataset, we compute the mean and standard deviation of the prediction MAEs from each fold of 10-fold cross validation in \textit{ML CV}. This MAE is different from our previous MAEs, denoted as \textit{ML average}, which uses the mean from cross validations to make a single prediction. Overall, \textit{ML average} outperforms a linear correction for all properties, indicating the generality of the customized correction of systematic errors. However, there is a relatively high variance between different folds of cross validation due to the small data size, especially for the extrapolation tasks. \textit{ML CV} performs the same or slightly worse than a linear correction for $D_{\mathrm{TFSI}}$, $D_{\mathrm{Poly}}$, and $D_{\mathrm{Poly}}^*$. A potential explanation is that a linear correction already performs reasonably well for these properties, demonstrated by the small MAEs of linear correction, while a more complicated multi-task model is prone to overfitting the noises in a small $\SI{50}{ns}$ dataset.

In Fig. \ref{fig:systematic-error}(c), we further study how the performance of our model would evolve with less $\SI{50}{ns}$ data, since these long MD simulations are expensive to run and cannot be easily parallelized. We find that the performance of the multi-task model decreases relatively slowly with less training data, and it still has some correction ability even with 13 CV data points, despite the large uncertainties due to the small data size. This obsevation shows the advantage of co-training a larger $\SI{5}{ns}$ dataset and a smaller $\SI{50}{ns}$ dataset -- it is much easier to learn a systematic correction than learn the property from scratch, and the co-training allows the transfer of graph representation learning from the $\SI{5}{ns}$ dataset to the $\SI{50}{ns}$ dataset. In contrast, the performance of a single-task model directly predicting $\SI{50}{ns}$ conductivity degrades much faster with less training data.

\begin{figure*}[tb]
  \centering
  \includegraphics[width=\linewidth]{exploration.pdf}
  \caption{Screening of polymer electrolytes. (a) Illustration of the screening workflow. (b) Distribution of the conductivities of top 50 polymers after each iteration. (c) Predictions of $\SI{50}{ns}$ and $\SI{5}{ns}$ conductivities for 6247 polymers in the search space. Green line denotes the top 50 conductivity from both predictions. (d) CPU hours that are actually used, required to screen the entire 6247 search space, and required to screen the 53362 candidate space.}
  \label{fig:exploration}
\end{figure*}

\textbf{Acceleration of the screening of polymers.} After demonstrating the performance of the multi-task model on reducing both random and systematic errors, we employ this model to perform an extensive screening of polymer electrolytes in the polymer space defined in Fig. \ref{fig:illustrative}(a). The goal of the screening is to search for polymers with the highest conductivity. As shown in Fig. \ref{fig:exploration}(a), we obtain 53362 polymer candidates using polymerization criteria from the ZINC chemical database \cite{irwin2005zinc}. To reduce the average computational cost, we limit our search space to only include polymers with monomer molecular weight less than 200, resulting in 6247 polymers.

We first use $\SI{5}{ns}$ MD simulations and a single-task GCN to explore polymers in the search space. To reduce computational cost, we only simulate each polymer once and employ GCN to reduce the random errors in the simulation. We perform 300 simulations in each iteration, 150 on randomly sampled polymers and 150 on best polymers predicted by GCN, which balances the exploration and exploitation. As shown in Fig. \ref{fig:exploration}(b), the conductivities of the top 50 polymers gradually increase as more polymers are explored with the iterative approach. But after 900 simulations, the average conductivity only increases slightly, indicating that we have achieved the best polymers in the 6247 search space based on $\SI{5}{ns}$ simulations.

Due to the systematic differences between $\SI{5}{ns}$ and $\SI{50}{ns}$ simulations, we randomly sample 120 polymers from those 900 polymers (876 successful simulations) and perform additional $\SI{50}{ns}$ MD, in which 117 are successful. These data allow us to correct the systematic errors in $\SI{5}{ns}$ simulation using the multi-task model. In Fig. \ref{fig:exploration}(c), we use the multi-task model to predict the $\SI{50}{ns}$ and $\SI{5}{ns}$ conductivities of all 6247 polymers in the search space. As a result of the customized correction, the ordering of conductivity changes from $\SI{5}{ns}$ to $\SI{50}{ns}$ predictions. The Spearman's rank correlation coefficient between these two predictions is 0.852, indicating that the ordering change is small but significant. For the top 50 polymers from $\SI{5}{ns}$ predictions, only 37 remain in the top 50 based on $\SI{50}{ns}$ predictions. This ordering change shows that the correction of systematic errors help us to identify some polymers that might be disregarded if only 5 ns simulations are performed.

To estimate the amount of acceleration we achieve, we compare the actual CPU hours used to the CPU hours that would be required if we performed one $\SI{50}{ns}$ MD simulations for each polymer. These simulations are run on NERSC Cori Haswell Compute Nodes and the CPU hours are estimated by averaging 100 simulations. In total, we use approximately 394,000 CPU hours for the MD simulations, with 33.2\% for sampling and relaxing amorphous structure, 28.6\% for $\SI{5}{ns}$ MD, and 38.2\% for $\SI{50}{ns}$ MD. The total cost only accounts for around 4.4\% and 0.51\% of the computation needed to simulate all the polymers from the 6247 search space and the 53362 candidate, respectively. Note that this conservative estimation assumes that only one $\SI{50}{ns}$ MD simulation is performed for each polymer for the brute-force screening. As shown in the previous section, our model has a true prediction error smaller than the random error from a $\SI{5}{ns}$ MD simulation. Although the random error from $\SI{50}{ns}$ simulation might be smaller, our model may have a larger acceleration due to the effect of random error reduction.

\textbf{Validation of the best candidates from the screening.} We employ the learned multi-task model to screen all 6247 polymers in the search space and 53362 polymers in the candidate space. In Fig. \ref{fig:validation}(a), we use \SI{50}{ns} MD to simulate 10 polymers out of the top 20 in the search space and 14 polymers out of the top 50 in the candidate space. These polymers are randomly selected from the top polymers using Butina clustering \cite{butina1999unsupervised,rdkit} to reduce their structural similarity, and only polymers which have not been seen in the \SI{50}{ns} dataset are selected. We observe a MAE of \SI{0.120}{log_{10}(S/cm)} and \SI{0.093}{log_{10}(S/cm)} for the predictions in search space and candidate space, respectively, which are between the interpolation and extrapolation errors in Fig. \ref{fig:systematic-error} and Table \ref{tab:comparison}. It shows that the extrapolation to the candidate space is easier than our hypothetical extrapolation test in Fig. \ref{fig:systematic-error}(b), yet a similar underestimation of conductivity is observed in the extrapolation. The larger errors for the top polymers in the search space might be explained by a combination of extrapolation errors and random errors in 50 ns MD simulations. We summarize the structure of the top polymers in Table \ref{tab:SI-interpolation} and \ref{tab:SI-extrapolation}, and most of them have PEO-like substructures which might explain their relative high conductivity.

In Fig. \ref{fig:validation}(b), we further validate the prediction of the model by gathering experimental conductivities for 31 different polymers from literature which are measured at the same salt concentration and temperature as our simulations \cite{pesko2016universal,zheng2018optimizing,tominaga2017ion,meabe2017polycondensation,mindemark2016allyl,fonseca2006development,itoh2018thermal,pehlivan2010ionic,he2017carbonate,pesko2016effect,doeff2000transport,silva2006novel}, and the results are also summarized in Table \ref{tab:SI-experiment}. Note that some polymers, like polyethylene oxide (PEO), do not follow the same structure pattern as our polymers. Nevertheless, the model still gives reasonable prediction on these out-of-distribution polymers because there are many PEO-like polymers in the training data. The largest errors come from the polymers with experimental conductivity less than $10^{-5}$ S/cm. In general it is difficult to simulate the conductivity of polymers with such low conductivity due to the long MD simulation time needed for convergence. In Fig. \ref{fig:SI-md-comparison}, we observe a much smaller prediction error with respect to 50 ns MD simulated conductivities for these polymers, indicating that the error with respect to the experiments are likely caused by the limited simulation time in MD. Other than the difficulty of simulating low-conductivity polymers, possible causes of the error also include the inaccuracy of the force fields, the finite length of the polymer chain, the finite size of the simulation box, etc. For the top polymers like PEO, we observe an underestimation of conductivity because the model cannot extrapolate to these polymers that are significantly different from the training data.

\begin{figure*}[tb]
  \centering
  \includegraphics[width=0.67\linewidth]{validation.pdf}
  \caption{Validation of the predicted polymers. (a) Validation of the best candidates from the search space (blue) and the candidate space (green). (b) Validation of the model prediction with out-of-distribution literature data.}
  \label{fig:validation}
\end{figure*}

\section*{Discussion}

We have performed a large scale computational screening of polymer electrolytes by learning to reduce random and systematic errors from molecular dynamics simulation with a multi-task learning framework. Our screening shows that PEO-like structure is the optimum structure for a broad class of carbonyl-based polymers. This finding is not surprising given that PEO has been one of the best polymer electrolyte since its discovery in 1973 \cite{fenton1973complexes}. However, with the large space of polymers screened in this work, our results reconfirm the advantage of PEO-like structures with respect to a broad range of chemical structures. Several potential directions remain open for discovering polymer electrolytes better than PEO. The first is to search for polymer electrolytes that achieve optimum conductivity at very high salt concentrations. Conductivity generally increases with increased salt concentration, but ion clustering and decreased diffusivity will reduce conductivity at high concentrations. \cite{hallinan2013polymer} Our screening keeps a constant concentration of \SI{1.5}{mol/kg} LiTFSI for different polymers, but some polycarbonate electrolytes show advantage at an extremely high salt concentrations \cite{tominaga2014fast,tominaga2015effect}. The second is to explore polymer chemistry beyond this study. Due to the limitations of the Monte Carlo procedure used to generate initial configurations, our simulations do not include polymers with aromatic rings. Recent studies propose the potential of polymers with high fragility and aromatic rings as polymer electrolytes due to the decoupling of ionic conductivity from structural relaxation \cite{agapov2011decoupling}. Backbones containing different lewis acidic heteroatoms or non-carbonyl based motifs could also lead to better polymer electrolytes. \cite{savoie2017enhancing}

The large scale screening is possible because we significantly reduce the computational cost of individual simulations by learning from imperfect data with the multi-task learning framework. The ability of neural networks to learn from noisy data is extensively studied in machine learning \cite{rolnick2017deep,arpit2017closer,han2018co} and has recently been applied to reduce the signal-to-noise ratio of band-excitation piezoresponse force microscopy \cite{borodinov2019deep} in materials science. Despite the wide use of graph neural networks in material discovery \cite{ahmad2018machine,back2019toward,back2019convolutional}, the random errors in training data are less studied, possibly because previous studies focus on simpler materials of which the random errors are much smaller. We show that random errors can be effectively reduced by learning a graph neural network across different chemistry even when the random error for each simulation is significant. It provides a potentially generalizable approach to accelerate the screening of complex materials whose structures can only be sampled from a distribution, e.g. amorphous polymers, surface defects, etc., because only one, instead of several, simulation needs to be performed for each material by adopting our approach.

The systematic error reduction demonstrated in this work is closely related to the transfer learning studies that aim to combine data from different sources \cite{yamada2019predicting,cubuk2019screening,smith2019approaching,zhu2020charting}. Our unique contribution in this work is to demonstrate the value of short, unconverged MD simulations in the context of material screening. We find that the systematic error between the \SI{5}{ns} and \SI{50}{ns} simulated transport properties can be corrected with a small amount of \SI{50}{ns} simulations, which can potentially be generalized to other types of materials, properties, and simulation methods. In summary, we hope that the random and systematic error reductions observed in this work could highlight the value of imperfect, cheaper simulations for material screening that might previously be overlooked. A broader class of complex materials could be screened with a similar approach if a cheap, noisy, and biased simulation method can be identified.

\section*{Methods}

\textbf{Graph representation for polymers.} The polymers are represented by graphs based on their monomer structure. The node embeddings $\bm{v}_i$ and edge embeddings $\bm{u}_{ij}$ are initialized using atom and bond features described in Table \ref{tab:SI-atom-feature} and \ref{tab:SI-bond-feature}. An additional edge is added to connect two ends of the monomer, allowing the end atoms to know the local chemical environments. We find that this representation has a better performance than using dummy atoms to denote the monomer ends.

\textbf{Network architecture.} We employ a graph convolution function developed in ref. \cite{xie2018crystal} to learn the node embeddings in the graph. For each node $i$, we first concatenate the center node, neighbor, and edge embeddings from last iteration $\bm{z}_{(i,j)}^{(t-1)} = \bm{v}_i^{(t-1)} \mathbin\Vert \bm{v}_j^{(t-1)} \mathbin\Vert \bm{u}_{(i, j)}$, then perform graph convolution,
\begin{equation}
	\bm{v}_i^{(t)} = \bm{v}_i^{(t-1)} + \sum_{j \in Neigh(i)} \sigma(\bm{z}_{(i,j)}^{(t-1)} \bm{W}_{\mathrm{f}}^{(t-1)} + \bm{b}_{\mathrm{f}}^{(t-1)}) \cdot g(\bm{z}_{(i,j)}^{(t-1)} \bm{W}_{\mathrm{s}}^{(t-1)} + \bm{b}_{\mathrm{s}}^{(t-1)}),
\end{equation}
where $\bm{W}_{\mathrm{f}}^{(t-1)}$, $\bm{W}_{\mathrm{s}}^{(t-1)}$, $\bm{b}_{\mathrm{f}}^{(t-1)}$, $\bm{b}_{\mathrm{s}}^{(t-1)}$ are weights, $\sigma$ and $g$ are sigmoid and softplus functions, respectively. After learning the node embeddings, we use a global soft attention pooling developed in ref. \cite{li2015gated} to learn a graph embeding,
\begin{equation}
	\bm{v}_{\mathcal{G}} = \sum_i \mathrm{softmax}(h_{\mathrm{gate}}(\bm{v}_i)) \cdot h(\bm{v}_i),
\end{equation}
where $h_{\mathrm{gate}}: \mathbb{R}^F \rightarrow \mathbb{R}$ and $h: \mathbb{R}^F \rightarrow \mathbb{R}^F$ are two fully connected neural networks. The graph embedding $\bm{v}_{\mathcal{G}}$ is then used in Eq. \ref{eq:network-5ns} and Eq. \ref{eq:network-50ns} to predict polymer properties.

\textbf{Molecular dynamics simulations.} The molecular dynamics simulations are performed with the large atomic molecular massively parallel simulator (LAMMPS) \cite{plimpton1995fast}. The atomic interactions are described by the polymer consistent force-field (PCFF+) \cite{sun1994force,rigby1997computer}, which has been previously used for polymer electrolyte systems \cite{molinari2018effect,france2019effect,france2019correlations}. The charge distribution of \ce{TFSI^-} is adjusted following ref. \cite{monteiro2008transport}, using a charge scaling factor of 0.7, to better describe the ion-ion interactions. All partial charges are reported in Table \ref{tab:SI-partial-charges}. There are 50 \ce{Li^+} and \ce{TFSI^-} in the simulation box. Each polymer chain has 150 atoms in the backbone. The number of polymer chains is determined by fixing the molality of LiTFSI at \SI{1.5}{mol/kg}. The initial configurations are generated using a Monte Carlo algorithm, implemented in the MedeA simulation environment \cite{MedeA}. The \SI{5}{ns} long equilibration procedure is based on a scheme described in ref. \cite{molinari2018effect}. Once equilibrated, the system is then run in the canonical ensemble (nVT) at a temperature of \SI{353}{K}, using a rRESPA multi-timescale integrator \cite{tuckerman1992reversible} with an outer timestep of 2 fs for non-bonded interactions, and an inner timestep of 0.5 fs. The high throughput workflow is implemented using the FireWorks workflow system \cite{jain2015fireworks}.

\textbf{Calculation of transport properties.} The diffusivities of lithium and TFSI ions are calculated using the mean squared displacement (MSD) of the corresponding particles,
\begin{equation}
	D = \frac{\left<[\bm{x}_i(t) - \bm{x}_i(0)]^2\right>}{6t},
\end{equation}
where $\bm{x}$ is the position of the particle, $t$ is the simulation time, and $\left<\cdot\right>$ denotes an ensemble average over the particles. The diffusivity of the polymer is calculated by averaging the diffusivities of \ce{O}, \ce{N}, and \ce{S} atoms in the polymer chains. The conductivity of the entire polymer electrolyte is calculated using the cluster Nernst-Einstein approach developed in ref. \cite{france2019correlations}. This method takes into account ion-ion interactions in the form of aggregation of ion clusters,
\begin{equation}
	\sigma = \frac{e^2}{Vk_BT} \sum_{i=0}^{N_+}\sum_{j=0}^{N_-}z_{ij}^2 \alpha_{ij}D_{ij},
\end{equation}
where $\alpha_{ij}$ is population of the ion clusters containing $i$ cations and $j$ anions, $z_{ij}$, $D_{ij}$ are the charge and diffusivity of the cluster, $N_+$ and $N_-$ are the maximum number of cations and anions in the clusters, $e$ is the elementary charge, $k_B$ is the Boltzmann constant, and $V$ and $T$ are the volume and the temperature of the system. We use the $\mathrm{cNE}_0$ approximation that assumes $D_{ij}$ is equal to the average diffusivity of lithium ion if the cluster is positively charged, and TFSI ion if the cluster is negatively charged. \cite{france2019correlations}.

\textbf{Data availability.} The toy LogP dataset, the \SI{5}{ns} and \SI{50}{ns} MD datasets are available at URL. The experimentally measured conductivity from literature is available at supplementary Table \ref{tab:SI-experiment}. The MD trajectories will be made available in a future publication.

\textbf{Code availability.} The multi-task graph neural network is implemented with PyTorch \cite{paszke2019pytorch} and PyTorch Geometric \cite{Fey/Lenssen/2019}. The code is available at URL.

\section*{Supplementary Notes}

\textbf{Supplementary Note 1: estimate true prediction error from noisy data.}

We assume there exists a deterministic function $f$ that maps from the polymer structure $\mathcal{G}$ to its true target property. However, due to the random errors associated with the initial configuration in MD simulations, the simulated target property $t$ has a small random error $\epsilon$,
\begin{equation}
	t = f(\mathcal{G}) + \epsilon,
\end{equation}
where $\epsilon$ follows a normal distribution with zero bias $\mathcal{N}(0, \eta)$. Here, we assume that $\epsilon$ is not a function of $\mathcal{G}$, i.e. different polymers have the same random error independent of their structure. This assumption is approximately correct based on the differences in conductivity of the same polymer between two independent MD simulations in the log scale (Fig.~\ref{fig:SI-true-error}).

To estimate the true prediction error of our model, we write our graph neural network model as a deterministic function $g$ that predicts polymer property based on their structure $\mathcal{G}$,
\begin{equation}
	y = g(\mathcal{G}).
\end{equation}
Note that we use different labels for the predicted property $y$ and the MD simulated property $t$.

Under these assumptions, the mean squared error between ML predictions and MD simulated properties, i.e. apparent prediction error, is,
\begin{equation}
	\mathrm{MSE}(y, t) = \mathop{\mathbb{E}}_\mathcal{G}[(y - t)^2] = \mathop{\mathbb{E}}_\mathcal{G}[(y - f(\mathcal{G}) - \epsilon)^2] = \mathop{\mathbb{E}}_\mathcal{G}[(y - f(\mathcal{G}))^2] + \mathop{\mathbb{E}}_\mathcal{G}[\epsilon^2].
\end{equation}
Note that in the last step we use the fact that $\mathop{\mathbb{E}}_\mathcal{G}[\epsilon] = 0$.

The mean squared error between two independent MD simulations for the same polymer is,
\begin{equation}
	\mathrm{MSE}(t_1, t_2) = \mathop{\mathbb{E}}_\mathcal{G}[(t_1 - t_2)^2] = \mathop{\mathbb{E}}_\mathcal{G}[(\epsilon_1 - \epsilon_2)^2] = 2\mathop{\mathbb{E}}_\mathcal{G}[\epsilon^2].
\end{equation}

Therefore, the mean squared error between ML predictions and the true target property, i.e. true prediction error, is,
\begin{equation}
	\mathrm{MSE}(y, f(\mathcal{G})) = \mathop{\mathbb{E}}_\mathcal{G}[(y - f(\mathcal{G}))^2] = \mathrm{MSE}(y, t) - \mathrm{MSE}(t_1, t_2) / 2.
\end{equation}

Based on our predictions on 86 testing data, $\mathrm{MSE}(y, t) = 0.0209$ and $\mathrm{MSE}(t_1, t_2) = 0.0274$. Therefore, the true prediction error $\mathrm{MSE}(y, f(\mathcal{G})) = 0.0072$. In comparison, the random error $\eta^2 \approx \mathop{\mathbb{E}}_\mathcal{G}[\epsilon^2] = 0.0137$. Remember that random errors can be reduced by running multiple MD simulations on the same polymers and computing the mean of target properties. Since $\eta_n = \eta / \sqrt{n}$, we estimate our ML prediction accuracy is approximately the accuracy of running $1.9 \approx 2$ MD simulations for each polymer.

\clearpage
\section*{Supplementary Figures}

\begin{figure*}[h]
  \centering
  \includegraphics[width=0.5\linewidth]{SI_md_50ns_exp.pdf}
  \caption{Comparison between 50 ns MD simulated conductivity and experimental conductivity reported in literature at the same salt concentration and temperature.}
  \label{fig:SI-md-50ns-exp}
\end{figure*}

\begin{figure*}[h]
  \centering
  \includegraphics[width=0.5\linewidth]{SI_true_error.pdf}
  \caption{Differences in conductivity of the same polymer between two independent 5 ns molecular dynamics simulations.}
  \label{fig:SI-true-error}
\end{figure*}

\begin{figure*}[h]
  \centering
  \includegraphics[width=0.5\linewidth]{SI_md_comparison.pdf}
  \caption{Comparison between 50 ns MD and predicted conductivity for the polymers from literature.}
  \label{fig:SI-md-comparison}
\end{figure*}

\clearpage
\section*{Supplementary Tables}

\begin{table*}[h
\caption{Molecular structure, simulated conductivity, and predicted conductivity for the top polymers in the search space.\label{tab:SI-interpolation}}
\begin{ruledtabular}
{\small\renewcommand{\arraystretch}{.8}
\begin{tabular}{p{7cm}cc}
  SMILES  & 50 ns MD $\sigma$ & Predicted $\sigma$ \\
  \hline
CN(CCCO[Cu])CCOCCOC(=O)[Au] & -3.74 & -3.95\\

O=C([Au])NCCOCCOCCOCCO[Cu] & -3.76 & -3.73\\

O=C([Au])OCCOCCCCOCCO[Cu] & -3.85 & -3.90\\

O=C([Au])OCCCSCCOCCO[Cu] & -3.97 & -3.95\\

O=C([Au])OCCCOCCCOCCCO[Cu] & -4.00 & -3.94\\

C=CCN(CCO[Cu])CCOCCOC(=O)[Au] & -4.00 & -3.95\\

O=C([Au])NCCOCCOCCOCCN[Cu] & -4.00 & -3.82\\

O=C([Au])OCCNCCOCCO[Cu] & -4.02 & -3.85\\

NOCCNCC(COC(=O)[Au])O[Cu] & -4.10 & -3.93\\

CN(CCO[Cu])CCOCCOC(=O)[Au] & -4.19 & -3.92\\

\end{tabular}}
\end{ruledtabular}
\end{table*}

\begin{table*}[h
\caption{Molecular structure, simulated conductivity, and predicted conductivity for the top polymers in the candidate space.\label{tab:SI-extrapolation}}
\begin{ruledtabular}
{\small\renewcommand{\arraystretch}{.8}
\begin{tabular}{p{10cm}cc}
  SMILES  & 50 ns MD $\sigma$ & Predicted $\sigma$ \\
  \hline
O=C([Au])OCCOCCOCCOCCOCCOCCOCCOCCOCCOCCOC
COCCOCCOCCOCCOCCOCCOCCOCCO[Cu] & -3.25 & -3.42\\

O=C([Au])NCCOCCOCCOCCOCCOCCOCCOCCOCCOCCOC
COCCOCCOCCOCCOCCOCCOCCOCCOCCN[Cu] & -3.32 & -3.48\\

O=C([Au])NCCOCCOCCOCCOCCOCCOCCOCCOCCOCCOC
COCCO[Cu] & -3.32 & -3.49\\

O=C([Au])OCCOCCOCCOCCOCCOCCOCCOCCOCCOCCOC
COCCOCCOCCN[Cu] & -3.36 & -3.47\\

O=C([Au])OCCOCCOCCOCCOCCNCCOCCOCCOCCOCCO[Cu] & -3.39 & -3.50\\

CN(CCOCCOCCOCCO[Cu])CCOCCOCCOCCOC(=O)[Au] & -3.39 & -3.54\\

C=CCOCC(CNCCOCCOCCOCCOCCOCCO[Cu])OC(=O)[Au] & -3.54 & -3.67\\

C=CCCC(CNCCOCCOCCOCCOCCOCCOC(=O)[Au])O[Cu] & -3.56 & -3.70\\

C=CC(CNCCOCCOCCOCCOCCOCCO[Cu])OC(=O)[Au] & -3.59 & -3.63\\

CC(CNCCOCCOCCOCCOCCOCCOC(=O)[Au])O[Cu] & -3.62 & -3.63\\

CCC(CNCCOCCOCCOCCOCCOCCOC(=O)[Au])O[Cu] & -3.65 & -3.67\\

CCOCC(CNCCOCCOCCOCCOCCOCCO[Cu])OC(=O)[Au] & -3.65 & -3.68\\

COCC(CNCCOCCOCCOCCOCCOCCO[Cu])OC(=O)[Au] & -3.66 & -3.61\\

CC(C)OCC(CNCCOCCOCCOCCOCCOCCO[Cu])OC(=O)[Au] & -3.70 & -3.70\\

\end{tabular}}
\end{ruledtabular}
\end{table*}

\begin{table*}[h
\caption{Molecular structure, experiment conductivity, and predicted conductivity for polymers from literature.\label{tab:SI-experiment}}
\begin{ruledtabular}
{\small\renewcommand{\arraystretch}{.8}
\begin{tabular}{p{10cm}cc}
  SMILES  & 50 ns MD $\sigma$ & Predicted $\sigma$ \\
  \hline
[Cu]CCOCCOCCOCCOCCOCC[Au] & -2.90 & -3.57\\

[Cu]OCC[Au] & -2.96 & -3.49\\

[Cu]CCCCCCOCCOCCOCCOCCOCC[Au] & -3.03 & -3.65\\

[Cu]CCCCOCCOCCOCCOCC[Au] & -3.12 & -3.64\\

[Cu]CCCCOCCOCCOCCOCCOCC[Au] & -3.18 & -3.61\\

[Cu]CCOCCOCCOCCOCC[Au] & -3.20 & -3.59\\

[Cu]CCCCCCOCCOCCOCCOCC[Au] & -3.23 & -3.69\\

[Cu]OCOCCOCC[Au] & -3.32 & -3.52\\

[Cu]OC(CCCCC[Au])=O & -3.57 & -4.49\\

[Cu]OC(=O)OCCOCCOCC[Au] & -3.57 & -3.75\\

[Cu]OC(=O)OCCOCCOCCOCCOCCOCC[Au] & -3.59 & -3.60\\

[Cu]OC(=O)OCCCCCCCCCCCC[Au] & -3.65 & -4.30\\

[Cu]CC(C)O[Au] & -3.81 & -4.15\\

[Cu]OC(=O)OCCCCCCC[Au] & -3.81 & -4.34\\

[Cu]OC(=O)OCCOCCOCCOCC[Au] & -3.84 & -3.67\\

[Cu]OC(=O)OCCCCCCCCCC[Au] & -3.92 & -4.31\\

[Cu]OC(=O)OCCCCCC[Au] & -4.06 & -4.35\\

[Cu]OC(=O)OCCCCCCCCC[Au] & -4.06 & -4.32\\

[Cu]OC(=O)OCCOCC[Au] & -4.07 & -3.95\\

[Cu]C(=O)CCCC(=O)OC(C)CO[Au] & -4.10 & -4.78\\

[Cu]OC(=O)OCCCC[Au] & -4.20 & -4.40\\

[Cu]OC(=O)OCCCCC[Au] & -4.31 & -4.37\\

[Cu]OC(=O)OCCCCCCCC[Au] & -4.34 & -4.33\\

[Cu]NCC[Au] & -4.50 & -3.84\\

[Cu]C(=O)COCC(=O)OC(C)CO[Au] & -4.74 & -4.40\\

[Cu]OC(=O)OCC(CC)(COCC=C)C[Au] & -4.87 & -4.72\\

[Cu]OC(=O)OCCC[Au] & -5.04 & -4.46\\

[Cu]OC(=O)OC(C)C[Au] & -5.34 & -4.90\\

[Cu]OC(=O)OC(CC)C[Au] & -6.34 & -4.89\\

[Cu]OC(=O)OC(CCC)C[Au] & -6.35 & -4.94\\

[Cu]OC(=O)OCC(OC)(OC)C[Au] & -6.86 & -4.59\\

\end{tabular}}
\end{ruledtabular}
\end{table*}

\begin{table}[h
\caption{Atom features.\label{tab:SI-atom-feature}}
\begin{ruledtabular}
\begin{tabular}{ll}
  Feature & Description\\
  \hline
  Atom type & Atomic number of elements (one-hot)\\
  Degree & Atom degree (one-hot)\\
  Formal charge & Formal charge of atoms (one-hot)\\
  Number of hydrogen & Number of connected hydrogen atoms (one-hot)\\
  Hybridization & Hybridization type of the atomic orbitals (one-hot)\\
  Aromatic & Whether the atom belongs to an aromatic ring (binary)\\
  Ring & Whether the atom belongs to a ring (binary)\\

\end{tabular}
\end{ruledtabular}
\end{table}

\begin{table}[h
\caption{Bond features.\label{tab:SI-bond-feature}}
\begin{ruledtabular}
\begin{tabular}{ll}
  Feature & Description\\
  \hline
  Bond type & Type of the bond, e.g. single, double, aromatic (one-hot)\\
  Stereochemistry & Stereochemistry of the bond (one-hot)\\
  Conjugated & Whether the bond is conjugated (binary)\\

\end{tabular}
\end{ruledtabular}
\end{table}

\begin{table}[h
\caption{Partial charges of the ions.\label{tab:SI-partial-charges}}
\begin{tabular}{cc} \hline\hline
  Species & Charge ($e$)\\
  \hline
  Li & $+0.7000$ \\
  S & $+0.3395$ \\
  C & $+0.2100$ \\
  F & $-0.0826$ \\
  O & $-0.2513$ \\
  N & $-0.2982$ \\
  \hline\hline
\end{tabular}
\end{table}

redpajama_set_name: RedPajamaArXiv