Detection of post-translational modification in a peptide with single-amino acid resolution using a graphyne nanopore: findings from molecular dynamics simulations

Abstract

The translocation of molecules through nanopores represents an established technology for molecular identification via their time-dependent ionic current signal. If single-molecule protein sequencing and the identification of post-translational modifications, such as phosphorylation, could be successfully realized with nanopores then this would represent a major breakthrough for applications in the field of life sciences. Toward this goal, we explored the ionic current sensitivity for non-phosphorylated (IEEEIYGEFD) and phosphorylated (IEEEIpYGEFD) forms of the amino acids sequences, using graphyne nanopore. Our study is based on molecular dynamics simulation with a classic force fields description of the system. We find that, from the ionic current, both forms of the peptide can be distinguished through their distinct ionic current traces. The results reveal that a graphyne membrane with an embedded nanopore represents a promising candidate for use in biosensors. Our discoveries on TYR modification and the relationship between ionic current and nanopore area support the viability of graphyne membranes as extremely sensitive nanosensors that can provide unique signatures for various biomolecules.

---

1 Introduction

Proteins form the building blocks of life and maintain the body's vital biological processes. It is therefore crucial for life sciences to have access to tools that enable the study of the protein contents of individual cells. However, with the currently available technology, it is a very difficult and time-consuming process to accurately identify individual proteins. Proteins are made of amino acids, which are chemically diverse, and therefore sequencing of proteins is a more difficult process compared to DNA sequencing. Conventional protein sequencing techniques such as Edman degradation¹ and mass spectrometry² may fail to disclose the complete sequence; however, biological and solid-state nanopores are emerging as viable alternatives for this task.

Solid-state nanopores are well-suited for the detection and manipulation of individual molecules because of their inherently narrow constriction regions. Additionally, they exhibit notable advantages, including chemical, thermal, and mechanical stability, as well as tunable size.³ Nanopore technology demonstrates significant potential for biological molecule detection, attributed to its high spatial resolution.^4–10 Nanopores in two-dimensional(2D) materials offer a distinct advantage owing to their molecular-scale thickness, facilitating the identification of individual amino acids or nucleotides during sequencing.¹¹ A variety of these materials has been suggested for the detection of individual molecules, including monolayer MoS₂,^12–15 bilayer MoS₂,¹⁶ graphene,^17,18 MXenes,¹⁹ and h-BN.^20,21 Graphyne is distinguished among 2D materials due to its unique properties, making it a promising candidate for various applications. This material exhibits significant semiconducting properties, characterized by an energy bandgap of 0.52 eV. The diverse carbon hybridizations in graphyne confer distinct electronic properties, rendering it a promising candidate for electronic device applications.^22,23

Research studies have been performed on DNA and protein sequencing techniques utilizing ion current detection via nanopores.^4,18,20,25 Previous experimental works^26–28 have utilized nanopore techniques to detect phosphorylation post-translational modifications, where changes in the ion current are measured to detect protein phosphorylation within long polypeptides that pass through a nanopore. A recent study²⁹ presents a detection mechanism for screening certain ions to stabilize and enhance anthrax nanopore electrical measurements. The authors demonstrate that the enhanced anthrax protective antigen protein serves as a suitable nanopore option for differentiating various categories of individual small compounds. Liang Wang et al. developed a unique nanopore single-molecule biosensor designed for label-free, ultrasensitive, stable, and user-friendly biomarker monitoring.³⁰ This methodology has promise for use in both fundamental research animal models and therapeutic settings. They have chosen the identification of Alzheimer's disease biomarkers as a case study to illustrate its effectiveness.

Molecular dynamics (MD) simulations enable the study of the microscopic kinetics of single-molecule translocation through nanopores, providing a powerful tool capable of atomically resolving the ionic transport of biomolecules driven by an electric field through a nanopore and the complex interaction between biomolecules and the nanopore walls.^7,20,31–33 In particular, Sathe et al. investigated the influence of pore area on detection efficiency by conducting MD simulations, varying the diameter of a graphene nanopore, and observed a dependence between resistance and pore area.³³ This finding aligns well with experimental results.^25,34 Similarly, Hu, Mao, and Ghosal³⁵ conducted MD simulations of a graphene sheet, varying the pore diameter within the range of 1.0 to 3.0 nm, and observed qualitative agreement with experimental results reported by Garaj et al.¹⁷ Their results demonstrate a relationship between conductance and nanopore diameter.

Here, we have used molecular dynamics simulations with classical force fields to investigate^24,36,37 IEEEIYGEFD and its post-translational modification in the form of phosphorylation (IEEEIpYGEFD) using 2D graphyne nanopores as a single molecule sensor. Firstly, the ionic current across the graphyne membrane was explored in two distinct nanopore sizes and shapes (triangular and hexagonal). We have also determined the ionic sensitivity as a function of the electric field strength for single peptide identification, in which the ionic current of the open nanopore was used as a reference. Our findings demonstrate that this graphyne nanopore setup enables the differentiation between two amino acid sequences that differ by only a single modification (i.e., phosphorylation at one amino acid site). More specifically, for the phosphorylated tyrosine (pY) modification, the sensitivity for detection was approximately 60%, while for tyrosine (Y), it was around 40%. In addition, we also analyzed the unoccupied area within the nanopore as a function of ionic current for both sequences, where a linear correlation between them is observed, enabling us to distinguish each sequence. Thus, our study provides an understanding of the application of solid-state nanopore technology for detecting alterations in amino acid sequences through ionic conductance.

2 Methodology

Molecular Dynamics (MD) simulations were performed using the GROMACS package.³⁸ Each simulation box consists of a single-layer graphyne membrane, an amino acid chain, water, and ions (Fig. 1a). The graphyne unit cell geometry was relaxed using first-principles calculations based on density functional theory (DFT). The membrane was generated by replicating the unit cell using atomsk code,³⁹ and a nanopore was created by removing selected atoms from the center of the graphyne supercell. Dangling bonds, created by the removal of carbon atoms, were saturated with hydrogen atoms. The dimensions of each box were adjusted to (48.0 × 47.5 × 82.0) Å, in x, y, and z, respectively. The membranes were fixed in xy plane at the center of the z-axis during all simulations.

Fig. 1 (a) Schematic illustration of the simulation system setup used for graphyne nanopore ionic current simulation. Water represented as semi-transparent blue volume; ions as yellow and purple sphere; the graphyne sheet in gray in the middle of the box; dangling bonds at the rim of the nanopore have been terminated by hydrogen atoms, shown in white. (b) Atomic model of peptide placed along the central axis of the nanopore. (c) Atomic model of the peptide simulated in this study, featuring the sequence IEEEIYGEFD and its phosphorylated form, IEEEIpYGEFD. This sequence belongs to the epidermal growth factor receptor (EGFR) protein.²⁴

The peptide chain was placed in the central region of the pore, aligning the center of mass of the amino acid currently under investigation with the z-axis position of the graphyne membrane. The simulation boxes were then filled with a 1.0 molar solution of KCl ions for all systems. To ensure that the amino acid currently under investigation remains at the center of the membrane for a sufficient amount of time to obtain a representative sample of the ionic current signal, one atom at each edge was fixed. The simulations were carried out for bias voltages of 1.0 V, 2.0 V, 3.0 V, and 4.0 V. The number of ions passing through the nanopore was recorded throughout the simulation and the ionic current was calculated. For each applied bias voltage we repeated the molecular dynamics simulations for different initial configurations, where the amino acid positions were changed, to ensure statistical sampling. The SPC/E⁴⁰ water model was utilized, and the Lincs algorithm was implemented to ensure the rigidity of the water molecules. The graphyne model employed the GAFF force field,⁴¹ with a 12.0 Å cutoff distance for Lennard-Jones interactions. Long-range electrostatic interactions were computed using the Particle–Mesh–Ewald method.⁴²

In each simulation, the steepest descent algorithm was applied over 50 000 steps to minimize energy and reduce repulsive contacts within the initial structures. Then, the system was equilibrated at the target temperature (300 K) for 10 nanoseconds (ns) using a constant pressure and temperature (NPT) ensemble. During this procedure, temperature, and pressure coupling were kept using a Velocity Rescaling thermostat^43,44 and semi-isotropic Parrinello–Rahman barostat⁴⁵ to the z direction, with coupling times set at 0.1 and 5.0 picoseconds (ps), respectively. Subsequently, a 5.0 ns simulation in the NVT ensemble was conducted employing the Nosé–Hoover thermostat, utilizing a coupling time of 0.5 ps. The production runs extended over 40 ns when an electric field was applied. To integrate the equations of motion, a leap-frog integrator with a femtosecond (fs) time step was employed. Atomic positions and energies were saved at regular intervals of 500 steps. Periodic boundary conditions were applied in all three orthogonal directions. Analyses of structural and dynamic properties were then performed using a combination of GROMACS³⁸ utilities and programs written by ourselves.

The ionic current through the graphyne nanopore as a function of time was calculated as the sum of all ions that cross the pore during the simulation time:^18,46

where δt was set to 10 ps, z_i(t) is the z-coordinate of ion i at time t, q_i is the charge of ion i, and N is the total number of ions.

In this paper, we named the amino acid sequence IEEEIYGEFD by assigning specific identifiers to each residue based on its position in the sequence. It consists of isoleucine (ILE1), followed by three glutamic acids (GLU2, GLU3, and GLU4), another isoleucine (ILE5), tyrosine (TYR6), glycine (GLY7), a fourth glutamic acid (GLU8), phenylalanine (PHE9), and aspartic acid (ASP10). When the tyrosine at position 6 is phosphorylated, it is denoted as phosphotyrosine (pTYR6), thus the phosphorylated form of the peptide is represented as IEEEIpYGEFD.

3 Results and discussion

To characterize the ionic conductance of graphyne nanopores, we performed classical MD simulations of the flow of ions through nanopores with two distinct diameters as a function of the applied external electric field, as shown in Fig. 1a and 2a. These simulations were conducted for pore A (B), with an average diameter of 1.26 nm (0.87 nm).

Fig. 2 (a) I–V curves of graphyne nanopores with pore diameters of 1.26 (pore A, shown in red color) and 0.87 nm (pore B, shown in blue color). The diameter was determined as the average distance between atoms on opposite sides of each pore; (b) Raw (gray) and LPF occupied-pore currents (in color) for simulations with various amino acids (as indicated by the legend in each graph in the upper right corner) positioned at the center of the nanopore, computed based on ion trajectories for IEEEIYGEFD peptide at a bias voltage of 4.0 V.

Analyzing the current–voltage curves for graphyne nanopore in Fig. 2a, we observed an Ohmic characteristic behavior. The conductance (G) values of 265.00 pS (73.15 pS) for nanopores with diameters of 1.26 nm (0.87 nm), were obtained. As expected, our results show that the conductance and nanopore size are directly related, which is consistent with previous works on similar systems.^17,25,33,35 We have in the following performed all simulations with pore A, with a larger diameter of 1.26 nm, since the smaller pore B, with its diameter of only 0.87 nm, might lead to difficulties in the translocation process of the peptide chain through the nanopore.

We now turn to the MD simulations of the amino acid sequence, both in its original form and when modified at one site through phosphorylation, when applying 1.0, 2.0, 3.0, and 4.0 V bias voltages. In these simulations, we aligned the z-coordinate center of mass of each amino acid under investigation with the nanopore plane. We have performed three simulations for each amino acid, slightly adjusting the center of mass of it to generate different samples and thereby improve the statistics. The initial positions of each amino acid, at the beginning of the simulations, are shown in Tables S1 and S2 (see SI). Using ion trajectories extracted from simulations, we calculated the ionic current using the eqn (1). The outcomes for various amino acids situated at the pore are shown in Fig. 2b for 4.0 V bias voltage (see Fig. S1–S4 for other bias voltage values). To minimize the noise in the calculated ionic current, we utilized a low-pass filter (LPF) with a cutoff frequency of 3.0 GHz.⁴⁷ Analyzing Fig. 2b, one can broadly discern that the ionic current average and its fluctuation are distinct between each amino acid, although it is not possible to quantify each one visually, and further analyses are therefore necessary.

The magnitude of the ionic currents was determined by averaging them over time and the samples. To quantify the resistance changes of each amino acid inside the nanopore the sensitivity was determined as follows:

where I_m is the ionic current of the target molecule and I_ref is the open pore current. The results for the voltage 4.0 V are presented in Fig. 3.

Fig. 3 Sensitivity of the ionic current for each amino acid for both IEEEIYGEFD and IEEEIpYGEFD peptides at a bias voltage of 3.0 and 4.0 V.

Analyzing the Fig. 3 for the IEEEIYGEFD and IEEEIpYGEFD peptides, we note that all amino acids present significant changes in their resistance, except for the GLU2 in the phosphorylated chain. We also observe that the ILE1, GLU4, GLY7, PHE9, and ASP10 have similar sensitivity considering both peptide chains (phosphorylated and non-phosphorylated). The amino acids GLU3, ILE5, TYR6, and GLU8 exhibit the most significant resistance change, comparing both peptide chains. Moreover, the sensitivity is higher for the non-phosphorylated chain than the phosphorylated one, except for the GLY7 amino (3.0 V). This fact is due to the change of the pore-area occupation shown in Fig. 6. These results indicate that the graphyne nanopore can discern differences between both sequences and that the modification of the amino acid TYR6 to pTYR6 impacts the ionic current not only associated with this amino acid but also with amino acids in the neighborhood up to the next-nearest neighbor. Furthermore, phosphotyrosine (pY) carries a negative phosphate group and consequently attracts K⁺. Residue-level radial distribution functions (RDF) show enhanced K⁺ proximity to pY relative to Y (Fig. 4a). A site-specific RDF to the phosphate oxygens (Fig. 4b) exhibits a pronounced first-shell peak, consistent with transient inner-sphere association of K⁺ near these oxygens. This increased K⁺ localization at pY is qualitatively consistent with the lower ionic current measured in pY-containing systems compared with Y. The values obtained with Δ_t = 50 ps are consistently greater than those with Δ_t = 10 ps; however, the current ratio between the molecule and open pore remains with a small deviation. This suggests that the selection of Δ_t does not affect the study's conclusions regarding the sensitivity and average current.

Fig. 4 (a) Radial distribution functions g(r) between ions and residue side chains: K–Y (blue), Cl–Y (green), K–pY (red), and Cl–pY (orange). (b) Site-specific g(r) between K⁺ and the phosphate oxygens of pY (label “O(PTYR)–K”).

Fig. 5 shows the histogram of the distributions of LPF ionic current for TYR6 and pTYR6 at different bias voltages. It is important to emphasize that our study compares two amino acid sequences in similar chemical environments. The unique distinction between them is the phosphorylation of the amino acid TYR6 to pTYR6. Despite this small difference, the LPF distributions present clearly visible deviations. One observes that for a small bias voltage, the distribution peaks at zero LPF ionic current. Increasing the voltage leads to the LPF current distributions shifting towards higher values. At 1.0 V bias voltage (Fig. 5a), the ionic current distributions for TYR6 and pTYR6 exhibit similar behavior, making it difficult to differentiate between them based on the ionic current. Furthermore, as noted above, it is evident that at this voltage bias, the distributions peak around zero. This outcome is attributed to the low voltage, which leads to only a few ions occasionally crossing the pore. For a bias voltage of 2.0 V (Fig. 5b), the distribution also peaks around zero, but with positive skewness⁴⁸ for both amino acids; however, this is somewhat more pronounced for pTYR6 than TYR6. Thus, for LPF ionic current higher than 200 pA the pTYR6 shows more events than TYR6.

Fig. 5 Histograms comparing LPF occupied-pore ionic currents for both peptides when tyrosine (Y) and phosphotyrosine (pY) are positioned at the center of the pore, for bias voltages of 1.0 V (a), 2.0 V (b), 3.0 V (c), and 4.0 V (d).

Conversely, at voltage biases of 3.0 and 4.0 V (Fig. 5c and d), a noticeable shift for higher values in the distributions is observed, indicating that the nanopore manifests distinct ionic currents for TYR6 and pTYR6. The events around zero decrease as expected for both considered voltages and the distributions show two soft peaks. Additionally, for the pTYR6 curves, zero values become less predominant, and the distributions are slightly different, allowing for possible identification. These results corroborate with the sensitivity discussed above.

To better understand the correlation between accessible pore area and ionic current, we have adapted the methodologies introduced by Heiranian, Farimani, and Aluru⁴⁹ and Liang et al.⁷ The SI contains more details about this methodology. It is important to emphasize that in our simulation, ILE1 and ASP10 amino acids are at the end of the chains, and the position of one atom of each one is constrained.

Fig. 6 shows the ionic current as a function of the unoccupied nanopore area. For the IEEEIYGEFD peptide, we verify a linear correlation between conductivity and available nanopore area, in agreement with previous results.⁷ The IEEEIpYGEFD peptide presents a similar behavior, but the slope of the curve decreases, and the conductivity increases compared to the previous one. This result is in agreement with the data shown in Fig. 5, where the distribution is seen to shift to higher values due to the lower pore blockage.

Fig. 6 Ionic current plotted as a function of the unoccupied pore area for an applied bias voltage of 4.0 V in the system with IEEEIYGEFD (a) and IEEEIpYGEFD (b). Green points denote the amino acids from GLU2 till PHE9, while the red point denote the amino acids ILE1 and ASP10. The red points were kept fixed, and were not included in the fitting.

4 Conclusions

We have theoretically studied the peptide IEEEIYGEFD and its phosphorylated version (IEEEIpYGEFD) using a graphyne membrane with a nanopore embedded. We employed a classical molecular dynamic simulation at 300 K and determined the ionic current trace for four different bias voltages (1.0 V, 2.0 V, 3.0 V, and 4.0 V). In light of our findings, the sensitivity of each considered amino acid is affected by bias voltage magnitude, identifying different groups of amino acids. Moreover, the amino acids ILE5, TYR6, pTYR6, and GLU8 present different sensitivities, leading to individual identification, specifically the capacity to identify the post-translational modification of amino acid TYR6 in the form of phosphorylation to pTYR6. The ionic current distributions support this result. Another significant analysis regarding the two amino acid chains (IEEEIYGEFD and IEEEIpYGEFD) was undertaken by considering the unoccupied nanopore area as a function of ionic current. In summary, our theoretical results lead us to the conclusion that solid-state nanopore sensors based on graphyne show great potential for identifying modifications in individual amino acid molecules.

Conflicts of interest

There are no conflicts to declare.

Data availability

Data will be made available on request. Supplementary information (SI) is available. See DOI: https://doi.org/10.1039/d5nr03335g.

Acknowledgements

The authors acknowledge financial support from the Brazilian agencies CAPES (001), FAPERJ, CNPq, FAPES (TO – 1043/2022), INCT Materials Informatics and computational resources of Sci-Com/ES, LNCC (SCAFMat2) and CENAPAD-SP as well as resources provided by the National Academic Infrastructure for Supercomputing in Sweden (NAISS), partially funded by the Swedish Research Council through grant agreement no. 2022-06725.

References

1. Y. Y. Elsayed, T. Kuhl and D. Imhof, Anal. Chem., 2024, 96, 4057–4066.
2. J. W. Robertson, C. G. Rodrigues, V. M. Stanford, K. A. Rubinson, O. V. Krasilnikov and J. J. Kasianowicz, Proc. Natl. Acad. Sci., 2007, 104, 8207–8211.
3. A. Nehra, S. Ahlawat and K. P. Singh, Sens. Actuators, B, 2019, 284, 595–622.
4. S. J. Heerema and C. Dekker, Nat. Nanotechnol., 2016, 11, 127–136.
5. Y. Feng, Y. Zhang, C. Ying, D. Wang and C. Du, Genomics, Proteomics Bioinf., 2015, 13, 4–16.
6. N. Yang and X. Jiang, Carbon, 2017, 115, 293–311.
7. L. Liang, P. Cui, Q. Wang, T. Wu, H. Ågren and Y. Tu, RSC Adv., 2013, 3, 2445–2453.
8. M. Zhang, C. Tang, Z. Wang, S. Chen, D. Zhang, K. Li, K. Sun, C. Zhao, Y. Wang and M. Xu, et al., Nat. Methods, 2024, 21, 609–618.
9. Z.-Y. Zhang, H.-L. Cui, D.-P. Huang and D.-Q. Wang, Sens. Actuators, B, 2021, 349, 130792.
10. X.-L. Xing, W. Li, L.-R. Guo, K. Wang, Y.-z. Ma, Q. Zhao and L. Ji, TrAC, Trends Anal. Chem., 2024, 117863.
11. Y. Goto, R. Akahori, I. Yanagi and K.-i. Takeda, J. Hum. Genet., 2020, 65, 69–77.
12. A. Barati Farimani, M. Heiranian and N. R. Aluru, npj 2D Mater. Appl., 2018, 2, 14.
13. M. Graf, M. Lihter, D. Altus, S. Marion and A. Radenovic, Nano Lett., 2019, 19, 9075–9083.
14. G. Sivaraman, F. A. De Souza, R. G. Amorim, W. L. Scopel, M. Fyta and R. H. Scheicher, J. Phys. Chem. C, 2016, 120, 23389–23396.
15. A. Perez, R. G. Amorim, C. E. Villegas and A. R. Rocha, Phys. Chem. Chem. Phys., 2020, 22, 27053–27059.
16. P. Sen and M. Gupta, RSC Adv., 2021, 11, 6114–6123.
17. S. Garaj, W. Hubbard, A. Reina, J. Kong, D. Branton and J. A. Golovchenko, Nature, 2010, 467, 190–193.
18. D. Zhao, H. Chen, Y. Wang, B. Li, C. Duan, Z. Li and L. Li, Front. Chem. Sci. Eng., 2021, 15, 922–934.
19. J. Prasongkit, S. Jungthawan, R. G. Amorim and R. H. Scheicher, Nano Res., 2022, 15, 9843–9849.
20. Z. Gu, Y. Zhang, B. Luan and R. Zhou, Soft Matter, 2016, 12, 817–823.
21. J. S. Lee, J. P. Oviedo, Y. M. N. D. Y. Bandara, X. Peng, L. Xia, Q. Wang, K. Garcia, J. Wang, M. J. Kim and M. J. Kim, Electrophoresis, 2021, 42, 991–1002.
22. N. Narita, S. Nagai, S. Suzuki and K. Nakao, Phys. Rev. B:Condens. Matter Mater. Phys., 2000, 62, 11146–11151.
23. D. C. Rodrigues, R. G. Amorim, A. Latgé and P. Venezuela, Carbon, 2023, 212, 118087.
24. T. Ohshiro, M. Tsutsui, K. Yokota, M. Furuhashi, M. Taniguchi and T. Kawai, Nat. Nanotechnol., 2014, 9, 835–840.
25. G. F. Schneider, S. W. Kowalczyk, V. E. Calado, G. Pandraud, H. W. Zandbergen, L. M. K. Vandersypen and C. Dekker, Nano Lett., 2010, 10, 3163–3167.
26. L. Restrepo-Pérez, C. H. Wong, G. Maglia, C. Dekker and C. Joo, Nano Lett., 2019, 19, 7957–7964.
27. W.-H. Lan, H. He, H. Bayley and Y. Qing, J. Am. Chem. Soc., 2024, 146, 24265–24270.
28. I. C. Nova, J. Ritmejeris, H. Brinkerhoff, T. J. Koenig, J. H. Gundlach and C. Dekker, Nat. Biotechnol., 2024, 42, 710–714.
29. J. Li, Y. Wang, L. Wang, Y. Wang, Z. Zhang, S. Liu and L. Wang, ACS Mater. Lett., 2025, 7, 2476–2481.
30. Q. Liu, Y. Ouyang, Y. Wang, S. Zhou, Y. Zhan and L. Wang, Adv. Healthcare Mater., 2025, 14, 2405058.
31. D. B. Wells, M. Belkin, J. Comer and A. Aksimentiev, Nano Lett., 2012, 12, 4117–4123.
32. H. Qiu, W. Zhou and W. Guo, ACS Nano, 2021, 15, 18848–18864.
33. C. Sathe, X. Zou, J.-P. Leburton and K. Schulten, ACS Nano, 2011, 5, 8842–8851.
34. C. A. Merchant, K. Healy, M. Wanunu, V. Ray, N. Peterman, J. Bartel, M. D. Fischbein, K. Venta, Z. Luo, A. T. C. Johnson and M. Drndić, Nano Lett., 2010, 10, 2915–2921.
35. G. Hu, M. Mao and S. Ghosal, Nanotechnology, 2012, 23, 395501.
36. T. Hunter and J. A. Cooper, Annu. Rev. Biochem., 1985, 54, 897–930.
37. S. V. Sharma, D. W. Bell, J. Settleman and D. A. Haber, Nat. Rev. Cancer, 2007, 7, 169–181.
38. D. Van Der Spoel, E. Lindahl, B. Hess, G. Groenhof, A. E. Mark and H. J. Berendsen, J. Comput. Chem., 2005, 26, 1701–1718.
39. P. Hirel, Comput. Phys. Commun., 2015, 197, 212–219.
40. H.-J.-C. Berendsen, J.-R. Grigera and T. P. Straatsma, J. Phys. Chem., 1987, 91, 6269–6271.
41. J. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman and D. A. Case, J. Comput. Chem., 2004, 25, 1157–1174.
42. U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577–8593.
43. L.-V. Woodcock, Chem. Phys. Lett., 1971, 10, 257–261.
44. G. Bussi, D. Donadio and M. Parrinello, J. Chem. Phys., 2007, 126, 014101.
45. M. Parrinello and A. Rahman, J. Appl. Phys., 1981, 52, 7182–7190.
46. A. Aksimentiev, J. B. Heng, G. Timp and K. Schulten, Biophys. J., 2004, 87, 2086–2097.
47. S. Cardoch, N. Timneanu, C. Caleman and R. H. Scheicher, ACS Nanosci. Au, 2022, 2, 119–127.
48. R. A. Groeneveld, Am. Stat., 1991, 45, 97–102.
49. M. Heiranian, A. B. Farimani and N. R. Aluru, Nat. Commun., 2015, 6, 8616.

---