A Hamiltonian replica swapping dynamics method and system for self-assembly systems

Abstract

This invention discloses a Hamiltonian replica exchange kinetics method for self-assembled systems, comprising: S1, determining whether the Hamiltonian replica exchange kinetics method is the REST method or the REST2 method; S2, implementing different Hamiltonian replica exchange kinetic transformations for the self-assembled system according to whether the Hamiltonian replica exchange kinetics method is the REST method or the REST2 method. The corresponding system and application are also disclosed. By reducing solute-solute and solute-solvent interactions between different replicas at a certain ratio, samples with different micelle distributions of the self-assembled system can be simulated more effectively. By adjusting appropriate parameter factors, the changes in solute-solute and solute-solvent interactions are roughly matched, thereby ensuring sufficient phase overlap between different replicas and the reliability of the entire simulation process, and enhancing the sampling efficiency of colloidal self-assembled system simulation.

Description

Technical Field

[0001] This invention belongs to the field of molecular dynamics simulation technology for biological genes and colloidal self-assembly, and particularly relates to a Hamiltonian copy exchange dynamics method and system for self-assembly systems. Background Technology

[0002] Amphiphilic molecules specifically refer to a class of functional molecules that simultaneously possess a hydrophilic head group and a lipophilic tail group, such as surfactants and phospholipids. When dissolved in water at concentrations exceeding the CMC (Closing Molar Concentration), these molecules tend to spontaneously aggregate due to the hydrophobic interaction of the tail group, forming micelles of varying sizes. The micelle self-assembly process is a crucial factor determining solution properties, and therefore has always been an important research direction with wide applications, including protein folding, colloidal emulsification, drug delivery, and the detergent industry. Experimentally, NMR and SCR techniques can analyze the structure and size distribution of micelles, but still lack sufficient microscopic information. With the improvement of computer hardware performance and the increasing maturity of technologies such as GPU parallel processing in recent years, simulating the self-assembly process of micelle systems using molecular dynamics methods and calculating the free energy differences between different states to examine the properties of these special systems from a microscopic perspective has become an increasingly important approach.

[0003] However, due to the specific nature of the system, the CMC concentration of amphoteric molecules is generally around 10. -3 Below mol / L, these concentrations fall into the category of dilute solutions, and are therefore limited by the most significant problem faced in general dilute solution simulations—excessively large spatial and temporal scales. Assuming that a dilute surfactant solution system can be abstractly divided into two states—aggregated and dispersed—to accurately simulate the system's properties, the simulation process must collect sufficient samples from both states. However, due to the extremely low concentration of the solution, in actual simulations, if the system initially is in a dispersed state, surfactant molecules will find it difficult to come into contact with other surfactants within the limited simulation time, thus hindering their transition to an aggregated state, and vice versa. Essentially, this is because a large free energy barrier exists between the two states of the system; overcoming this barrier to transition from one state to another is a very difficult and low-probability event.

[0004] To address the challenge of simulating low-probability events, one approach is to reduce the system's degrees of freedom, such as by employing implicit water models or coarse-grained force fields, thereby increasing the simulation timescale. However, this only works for specific systems (with relatively large CMCs) and is insufficient for most amphoteric molecular systems. Another approach is to use enhanced sampling methods to improve simulation efficiency. This approach can be broadly categorized into two types: one is biased simulations such as umbrella sampling, metastable dynamics, and accelerated molecular dynamics, which accelerate sampling efficiency by introducing biased forces or potential energy into the simulation to push the system over energy barriers. However, these methods require defining reaction coordinates before applying biased forces or potential energy along those coordinates, and currently, no suitable reaction coordinates exist to accurately describe surfactant self-assembly. The other type of enhanced sampling method can be classified as generalized ensemble methods, most commonly replica exchange dynamics. Since these methods do not depend on the selection of reaction coordinates, they are a better choice for self-assembling systems. However, with increasing system size, the relative fluctuations of total system energy decrease drastically in traditional replica exchange kinetics methods, leading to insufficient energy overlap between replicas and reduced replica exchange efficiency. To address this issue, Bernerd et al. developed the REST and REST2 methods. The basic idea behind these two methods is to selectively raise the temperature of the solute (e.g., protein) in the system while maintaining the temperature of the majority of the solution constant. This allows for significantly faster solute sampling in solution while ensuring replica exchange efficiency. Essentially, REST and REST2 are variations of Hamiltonian replica exchange kinetics, but a comprehensive Hamiltonian replica exchange kinetics scheme for self-assembling systems has not yet been developed.

[0005] Therefore, there is indeed a need to propose better solutions for the aforementioned existing technologies. Summary of the Invention

[0006] To address the shortcomings of existing technologies, a Hamiltonian replica exchange kinetics method and system for self-assembled systems is provided. By reducing solute-solute and solute-solvent interactions among different replicas at a certain ratio, samples with different micelle distributions in the self-assembled system can be simulated more effectively. By adjusting appropriate parameter factors, the changes in solute-solute and solute-solvent interactions are roughly matched, thereby ensuring sufficient phase overlap between different replicas and the reliability of the entire simulation process, and enhancing the sampling efficiency of colloidal self-assembled system simulation.

[0007] This invention provides a Hamiltonian copy exchange dynamics method for self-assembled systems, comprising:

[0008] S1, determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method;

[0009] S2, based on the Hamiltonian replica exchange dynamics method, implements different Hamiltonian replica exchange dynamics transformations for the self-assembled system using either the REST method or the REST2 method.

[0010] Preferably, the Hamiltonian replica exchange dynamics transformation of the self-assembled system according to the Hamiltonian replica exchange dynamics method for the REST method includes:

[0011] S21, determine the temperature of the m-th replica as T. m The biased potential function in this copy is:

[0012]

[0013] Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy;

[0014] S22, determine the scaling factor λ in different replicas; and respectively for the interaction E pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the following equation (4) as the potential function under temperature T0:

[0015]

[0016] S23 transforms each replica in REST into a simulation at the same temperature T0, with the only difference being the biased potential function.

[0017] Preferably, λ is greater than or equal to 0.

[0018] Preferably, the

[0019] Preferably, the proportionality coefficient λ is a factor. Therefore, equation (4) is transformed into the following potential function:

[0020]

[0021] Therefore, the m-th replica can also be regarded as a system simulated by Equation (2) as the potential function under the condition of temperature T0.

[0022] Preferably, the Hamiltonian replica exchange dynamics transformation of the self-assembled system for the REST2 method is performed according to the Hamiltonian replica exchange dynamics method, including:

[0023] S21, determine the temperature of the m-th replica as T. m The biased potential function in this copy is:

[0024]

[0025] Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy;

[0026] S22, in different replicas, determine the scaling factor λ; respectively for the interaction E pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the following equation (5) as the potential function under temperature T0:

[0027]

[0028] S23 transforms each replica in REST2 into a simulation at the same temperature T0, with the only difference being the biased potential function.

[0029] Preferably, λ is greater than or equal to 0.

[0030] Preferably, f(λ) is

[0031] Preferably, the proportionality coefficient λ is a factor. Therefore, equation (5) is transformed into the following potential function:

[0032]

[0033] A second aspect of the present invention provides a Hamiltonian copy exchange dynamics system for self-assembly systems, comprising:

[0034] The method determination module (101) is used to determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method;

[0035] The dynamic transformation module (102) is used to implement different Hamiltonian replica exchange dynamic transformations of the self-assembled system according to the Hamiltonian replica exchange dynamics method for the REST method or the REST2 method.

[0036] A third aspect of the invention is to provide the application of the Hamiltonian copy exchange kinetics method for self-assembled systems in the dynamic simulation of colloidal systems.

[0037] The method, system, and application provided by this invention have the following beneficial technical effects:

[0038] This invention, based on Hamiltonian replica exchange dynamics, extends the REST and REST2 methods to a more general case and applies them to molecular dynamics simulations of colloidal self-assembly. Specifically, monomer-monomer and monomer-water interactions are reduced to varying degrees along different replicas, allowing for sufficient sampling between polymeric and dispersed states. By appropriately setting the scaling parameter f and coupling parameter λ, sufficient phase space overlap and a large exchange probability can be controlled among all adjacent replicas, ensuring the high efficiency of replica exchange dynamics. Specific practical examples demonstrate that the method of this invention has higher simulation efficiency for the dynamics simulation of colloidal systems. Firstly, this method can obtain convergence results consistent with ordinary molecular dynamics simulations, including the maximum stable size of the colloid and the colloidal size distribution at equilibrium, and the computation time required to obtain the convergence results using this method is less. Attached Figure Description

[0039] Figure 1 This is a schematic diagram of the Hamiltonian replica exchange dynamics method for self-assembly systems according to a preferred embodiment of the present invention;

[0040] Figure 2 For the coarse-grained model of SNS surfactant and water shown in the preferred embodiment of the present invention, CT and CM represent carbon chains, SU represents sulfate, W represents water, and SOD represents hydrated sodium ions.

[0041] Figure 3 The image shows a screenshot of the initial conformation of a simulated system according to a preferred embodiment of the present invention. Water molecules are not shown in the image, the chain-like structures are SNS surfactants, and the single spheres are sodium ions.

[0042] Figure 4 This is a schematic diagram illustrating the HREMD simulated micellar polymerization process according to a preferred embodiment of the present invention;

[0043] Figure 5 This is a schematic diagram illustrating the variation of the average acceptance probability with the scaling parameter f according to a preferred embodiment of the present invention;

[0044] Figure 6 This is a schematic diagram illustrating the replica exchange diffusion ratio according to a preferred embodiment of the present invention;

[0045] Figure 7 This is a schematic diagram illustrating the motion changes of the replica in Hamiltonian space during the first 40 ns of simulation according to a preferred embodiment of the present invention;

[0046] Figure 8 This is a schematic diagram illustrating the potential energy distribution of different replicas and the potential energy overlap between different replicas according to a preferred embodiment of the present invention.

[0047] Figure 9 This is a schematic diagram illustrating the monomer number distribution of the surfactant solution in different simulated versions according to a preferred embodiment of the present invention.

[0048] Figure 10 This is a schematic diagram illustrating the change in the number of micelles in the surfactant solution over simulation time in HREMD and BFMD simulations according to a preferred embodiment of the present invention.

[0049] Figure 11 The illustration shows the solution micelle size distribution obtained by simulation calculation using HREMD and BFMD respectively, according to a preferred embodiment of the present invention. The shaded area represents the calculation deviation. In order to more clearly compare the relative sizes, the deviation of the calculation results of the two methods is listed separately in the illustration.

[0050] Figure 12 (A) and Figure 12 (B) is an HREMD ( ) shown according to a preferred embodiment of the present invention. Figure 12 A) and BFMD Figure 12 B) A schematic diagram of the autocorrelation coefficient and correlation time of the number of monomers in the simulation;

[0051] Figure 13 In the HREMD (above) and BFMD (below) simulations shown in the preferred embodiment of the present invention, the parameter χ 2 (t) Schematic diagram showing the change of t over time;

[0052] Figure 14 This is a diagram illustrating the Hamiltonian copy exchange dynamics system architecture for a self-assembly system according to a preferred embodiment of the present invention. Detailed Implementation

[0053] The specific embodiments of the present invention will be described in further detail below with reference to the accompanying drawings and examples. The following examples are for illustrative purposes only and are not intended to limit the scope of the invention.

[0054] The terminology used in the embodiments is as follows:

[0055] CMC: Critical micelle concentration

[0056] REMD: Replica Swapping Dynamics

[0057] HREMD: Hamiltonian Replica Swapping Dynamics

[0058] TREMD: Temperature Replica Exchange Dynamics

[0059] BFMD: Classical Molecular Dynamics

[0060] REST: Solute Heating Replica Exchange Kinetics

[0061] In self-assembling solution systems, during the aggregation of monomer molecules to form micelles, monomer-monomer pair interactions increase, while monomer-water pair interactions decrease, and there is a certain proportional relationship between the two. Therefore, by selecting an appropriate scaling factor, we can simultaneously adjust the relative magnitudes of these two contributions, ensuring that the energy changes of the two parts cancel each other out as much as possible. After adjusting the biased potential functions in different replicas using the above method, theoretically, the total energy difference between different replicas can be minimized, increasing the overlap between replicas. Ideally, even if the system conformations and micelle size distributions of different replicas vary greatly, as long as there is sufficient energy overlap between replicas, the exchange efficiency between replicas and the sampling efficiency of the system in the entire phase space will be greatly improved.

[0062] Both REST and REST2 essentially accelerate system kinetics and thus sampling by reducing solute-solute and solute-solution interactions across different replicas. In the REST method, let the temperature of the m-th replica be T. m The biased potential function in this copy can be written in the following form:

[0063]

[0064] Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th replica; the entire potential function is multiplied by a factor. The following potential function is obtained:

[0065]

[0066] Therefore, the m-th replica can also be regarded as a system simulated at temperature T0 with Equation (2) as the potential function. Thus, each replica in REST can be transformed into a simulation at the same temperature, the only difference being the biased potential function.

[0067] In REST2, all replicas are simulated at a target temperature, and the biased potential function of the m-th replica is as follows:

[0068]

[0069] Comparing equations (2) and (3), both involve changing E by a certain proportion. pp and E pw Interact, and maintain E ww The only difference is that it applies to E. pw The percentages of change are different.

[0070] If we consider all temperature effects as modifications to the potential function, then REST and REST2 are essentially Hamiltonian replica swapping dynamics methods. Assume a coupling parameter... Substituting into equations (1) and (2), we can obtain the following for REST and REST2 respectively:

[0071]

[0072]

[0073] These two equations show that in any copy, E pp The interaction changes proportionally to λ, while E pw (X) The interaction changes proportionally to f(λ), where f(λ) is in REST. In REST2 it is In fact, more generally speaking, as long as E is guaranteed pp and E pw The mutual balance between them ensures sufficient exchange efficiency between replicas, and f(λ) can be an arbitrary function. That is, by adjusting f(λ) more flexibly, we may be able to obtain higher replica exchange efficiency than conventional REMD, or even REST and REST2, thus enabling us to use this method to perform molecular dynamics simulations of colloidal self-assembly systems.

[0074] Example 1

[0075] See Figure 1 This embodiment provides a Hamiltonian replica exchange dynamics method for self-assembled systems, including:

[0076] S1, determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method;

[0077] S2, based on the Hamiltonian replica exchange dynamics method, implements different Hamiltonian replica exchange dynamics transformations for the self-assembled system using either the REST method or the REST2 method.

[0078] As a preferred embodiment, the Hamiltonian replica exchange dynamics transformation of the self-assembled system is implemented for the REST method according to the Hamiltonian replica exchange dynamics method, including:

[0079] S21, determine the temperature of the m-th replica as T. m The biased potential function in this copy is:

[0080]

[0081] Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy;

[0082] S22, determine the scaling factor λ in different replicas; and respectively for the interaction E pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the following equation (4) as the potential function under temperature T0:

[0083]

[0084] S23 transforms each replica in REST into a simulation at the same temperature T0, with the only difference being the biased potential function.

[0085] In a preferred embodiment, λ is greater than or equal to 0, which means that in different copies, surfactant-surfactant interactions are weakened by λ, while surfactant-water interactions are weakened by f(λ). That is, solute-solute and solute-solution interactions are weakened, thereby accelerating the kinetics of the system and thus accelerating the sampling of the system.

[0086] As a preferred embodiment, the

[0087] Of course, those skilled in the art will know that f(λ) can be any function. That is, by adjusting f(λ) more flexibly, it is possible to obtain a higher copy exchange efficiency than conventional REMD, or even REST and REST2, thereby enabling the use of this method to perform molecular dynamics simulations of colloidal self-assembly systems.

[0088] In a preferred embodiment, the proportionality coefficient λ is a factor. Therefore, equation (4) is transformed into the following potential function:

[0089]

[0090] Therefore, the m-th replica can also be regarded as a system simulated by Equation (2) as the potential function under the condition of temperature T0.

[0091] As a preferred embodiment, the Hamiltonian replica exchange dynamics transformation of the self-assembled system for the REST2 method is performed according to the Hamiltonian replica exchange dynamics method, including:

[0092] S21, determine the temperature of the m-th replica as T. m The biased potential function in this copy is:

[0093]

[0094] Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy;

[0095] S22, in different replicas, determine the scaling factor λ; respectively for the interaction W pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the potential function of equation (5) at temperature T0:

[0096]

[0097] S23 transforms each replica in REST2 into a simulation at the same temperature T0, with the only difference being the biased potential function.

[0098] In a preferred embodiment, λ is greater than or equal to 0, which means that in different copies, surfactant-surfactant interactions are weakened by λ, while surfactant-water interactions are weakened by f(λ). That is, solute-solute and solute-solution interactions are weakened, thereby accelerating the kinetics of the system and thus accelerating the sampling of the system.

[0099] In a preferred embodiment, f(λ) is

[0100] Of course, those skilled in the art will know that f(λ) can be any function. That is, by adjusting f(λ) more flexibly, it is possible to obtain a higher copy exchange efficiency than conventional REMD, or even REST and REST2, thereby enabling the use of this method to perform molecular dynamics simulations of colloidal self-assembly systems.

[0101] In a preferred embodiment, the proportionality coefficient λ is a factor. Therefore, equation (5) is transformed into the following potential function:

[0102]

[0103] Specific application examples:

[0104] I. Model Selection and Parameter Settings

[0105] Based on the method of this embodiment, a common anionic surfactant, sodium nonyl sulfate (SNS), was selected as the research object to test the effectiveness of the method.

[0106] Figure 2 A coarse-grained force field model of SNS molecules and water is presented, where CT and CM represent carbon chains, SU represents sulfate, W represents water, and SOD represents hydrated sodium ions. The initial system size is 10 × 10 × 10 nm. 3 The cubic box contains 11,068 coarse-grained water spheres and 400 coarse-grained SNS molecules. The concentration of the system is approximately 0.664 mol / L, higher than the CMC concentration of SNS. Initially, the surfactant molecules are randomly dispersed in the aqueous solution, such as... Figure 3 As shown, Figure 3 The image shows a screenshot of the initial conformation of the simulation system. Water molecules are not shown in the image, the chain-like structures are SNS surfactants, and the single spheres are sodium ions.

[0107] For HREMD simulations, the interaction E is performed in different replicas. pp and E pwMultiply by the proportionality constants 1-λ and 1-f(λ). Generally, λ is greater than or equal to 0, which means that in different copies, the surfactant-surfactant interaction is weakened by λ, while the surfactant-water interaction is weakened by f(λ) (e.g., Figure 4 As shown, Figure 4 (This is a schematic diagram of the micellar polymerization process simulated by HREMD). In addition, the reduction in surfactant-water interactions has an additional effect: it allows surfactant molecules to move faster in solution, accelerating the system's kinetics and contributing to improved sampling efficiency. In this example, eight replicas were established, with λ values ​​for each replica pre-tested at (0, 0.04, 0.0792, 0.1176, 0.1552, 0.192, 0.228, 0.2632) to ensure consistent exchange probabilities among different replicas. The parameter f was adjusted from 0.2 to 0.5 at 0.05 intervals, and the optimal value was found by comparing the average acceptance probabilities. Once parameter f was determined, all HREMD simulations used this parameter setting.

[0108] As a reference for testing simulation efficiency, a long-duration classical molecular dynamics simulation (Brute Force MD, BFMD) was also performed on the system. The efficiency of both HREMD and BFMD simulations was mainly evaluated by statistically analyzing the standard deviation of 10 parallel simulations. For a single sample, HREMD took 350 ns to simulate, while BFMD took 2800 ns, exactly eight times the simulation time of HREMD.

[0109] II. Results Display

[0110] (I) HREMD Simulation Efficiency Analysis

[0111] First, we need to find a suitable adjustment coefficient f so that, given a value of λ, the probability of swapping between replicas is maximized. Figure 5 The diagram illustrates how the average acceptance probability varies with the adjustment coefficient f. It shows that the average acceptance probability in the HREMD simulation peaks at f = 0.4; further increases or decreases in f lead to a decrease in the acceptance probability. This is because the scaling parameter f is essentially used to balance the opposing contributions of monomer-monomer interactions and monomer-water interactions. Therefore, if a suitable scaling parameter f is chosen so that these two contributions are essentially equal, the difference in total system energy between different replicas can be controlled within a very small range. This ensures maximum energy overlap between adjacent replicas and a higher replica exchange acceptance probability.

[0112] The replica diffusion ratio can be used to examine replica exchange efficiency. If we define the replica as being at a "low temperature" in the original potential function, the replica's "temperature" continuously increases as monomer-monomer and monomer-water interactions decrease. Therefore, the diffusion ratio can be defined as the ratio of the number of conformations just obtained from the "lowest temperature" to the total number of conformations just obtained from the "lowest temperature" and the "highest temperature". Matthias Troyer et al. have shown that in TREMD, when the temperature diffusion ratio is linearly related to the temperature index, the replica diffuses completely freely in the temperature space. Similarly, in HREMD, when the replica diffusion ratio is linearly related to the Hamiltonian "temperature" index, the replica diffuses most efficiently in the Hamiltonian space. Figure 6 The simulation results show the movement of replicas in Hamiltonian space during the first 40 ns of the simulation. The solid line describes the replica diffusion ratio of an arbitrarily selected HREMD simulation, while the dashed line represents the ideal straight line. Observation reveals that the difference between the results of this embodiment and the ideal state is very small, indicating that replicas can diffuse freely in the HREMD simulation without a significant bottleneck effect. Using the same sample as an example, the variation of a certain conformation's exchange between replicas over simulation time is also shown (e.g., Figure 7 (As shown). This only lists the case for the first 40ns of the simulation. It can be seen that the conformation shuttles back and forth between different replicas many times, which also proves the high replica switching efficiency from another perspective.

[0113] Figure 8 The system energy distribution across different replicas was shown. Observations revealed that the total energy distribution across different replicas generally conforms to an equally spaced Gaussian distribution, and the overlap between adjacent replicas is also largely consistent. This further confirms the equality of replica exchange probabilities in the HREMD simulation. Furthermore, using the number of monomers in the solution as an observation value also allows us to depict the distribution of that observation value across different replicas, such as... Figure 9 As shown, the distribution of monomer numbers in different replicas also exhibits a Gaussian distribution. The peak value of the Gaussian distribution monotonically increases along the replica sequence, indicating that the original potential function tends to aggregate the system, while the potential function after the interaction tends to disperse the system. Closer observation also reveals that as the system tends to disperse, the width of its distribution also increases, indicating that the system samples a wider distribution in phase space. Figure 8 This refers to the potential energy distribution of different replicas and the potential energy overlap between different replicas. Figure 9 The distribution of monomer number in surfactant solutions in different replica simulations.

[0114] (II) Comparison of HREMD and BFMD

[0115] Figure 10This describes the change in the total number of polymers in the system over simulation time, specifically the change in the number of micelles in the surfactant solution during HREMD and BFMD simulations. HREMD dots represent the HREMD method, with the simulation time displayed on the upper X-axis; BFMD diamonds represent the BFMD method, with the simulation time displayed on the lower X-axis, and its scale is exactly 8 times that of HREMD. Observations show that both methods initially decrease rapidly, then slowly reach an equilibrium value. This indicates that the polymerization of the surfactant in the system mainly consists of two stages: first, the polymerization of monomers into smaller oligomers is a rapid process; then, the oligomers further polymerize into micelles, and finally, the small micelles merge into larger micelles—a relatively slow process. Importantly, both methods eventually converge to similar regions, with an average polymer number between 11 and 12, represented by gray shading. This demonstrates that HREMD simulations can achieve similar convergence results to long-duration BFMD simulations, without requiring a longer total computation time. It is also worth noting that the standard statistical deviation of each data point in the HREMD results is smaller than that in the BFMD simulation, indicating that under the same computing time, HREMD can obtain more accurate results than BFMD. Here, 100 ns and 800 ns are selected as the equilibrium times for HREMD and BFMD, respectively, and all samples after this time are used to calculate the equilibrium properties of the system.

[0116] Figure 11 The micelle size distributions of the system simulated by BFMD and HREMD methods are presented respectively. Overall, both methods exhibit characteristics similar to experimentally measured micelle size distributions. A small number of surfactant molecules exist in solution as monomers or oligomers (dimeric, trimeric, etc.), while larger polymers (n in the range of 5 to 20) are almost absent. The main surfactant components are distributed in micelles ranging from 20 to 50. Comparing the two distributions, the peak positions are relatively consistent, both around 44. The peak height of the BFMD distribution is higher, while the width of the HREMD distribution is wider, indicating that HREMD collects more samples of different conformations. Of course, considering statistical bias (shown by shading in the figure), the micelle size distributions calculated by the two methods are considered to have a certain similarity, which also confirms the previous conclusion that HREMD simulations can obtain results similar to long-term BFMD simulations. The inset in the figure directly compares the statistical bias of the two methods. Clearly, HREMD has a smaller bias than BFMD across all ranges, especially in the range of large-sized micelle distributions.

[0117] (III) HREMD Sampling Efficiency Analysis

[0118] The amount of surfactant monomers in a system is an important characteristic of the system's conformation, and its changes also reflect the system's kinetic behavior. To verify the sampling efficiency of HREMD, the autocorrelation coefficient of the number of monomers in the solution was calculated, such as... Figure 12 As shown in the figure, the autocorrelation coefficient of HREMD decays significantly faster than that of BFMD. For quantitative comparison, integrating the autocorrelation coefficient to infinity yields the correlation time, as shown by the rising curve in the figure. The correlation times calculated by BFMD and HREMD are 1.2 ns and 0.13 ns, respectively. Considering that HREMD has a replica count of 8 (0.13 × 8 = 1.04, still less than 1.2), the shorter correlation time of HREMD compared to BFMD indicates that in simulations with equal computational time, HREMD can generate more independent samples than BFMD, further demonstrating that HREMD has a higher sampling efficiency than BFMD.

[0119] Finally, a parameter X is defined. 2 (t) is used to evaluate the convergence of simulation results for different samples.

[0120] χ 2 (t) describes how the difference in average micelle size distributions obtained from multiple parallel simulations changes over simulation time, and can be calculated using the following formula:

[0121]

[0122] Where m and n represent the number of parallel simulations and the specific micelle size, respectively, H i，j Specifically, this refers to the proportion of micelles of size j obtained statistically in the l-th simulation. Theoretically, as the simulation time increases, the system gradually reaches ergodicity, X 2 (t) should gradually approach 0. Figure 13 The variation of this parameter over time is described. Similarly, for ease of comparison, X-axis coordinates with different scales (8 times different) are set for BFMD (upper part) and HREMD (lower part). Comparing the changes of this parameter in the simulations of the two methods, it is found that the HREMD method corresponds to X... 2 (t) is always much smaller than the χ corresponding to the BFMD method. 2 (t). If 0.00001 is used as the convergence criterion, then HREMD and BFMD require 110ns and 3580ns respectively to converge. Considering the total number of replicas is 8, the HREMD method is still nearly 4 times faster than BFMD, thus further proving the high sampling efficiency of HREMD.

[0123] Example 2

[0124] See Figure 14A Hamiltonian copy swapping dynamics system for self-assembly systems includes:

[0125] Method determination module 101 is used to determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method;

[0126] The dynamic transformation module 102 is used to implement different Hamiltonian replica exchange dynamic transformations of the self-assembled system according to the Hamiltonian replica exchange dynamics method for either the REST method or the REST2 method.

[0127] Example 3

[0128] Application of Hamiltonian copy exchange dynamics method for self-assembled systems in the dynamic simulation of colloidal systems.

[0129] Although preferred embodiments of the invention have been described, those skilled in the art, upon learning the basic inventive concept, can make other changes and modifications to these embodiments. Therefore, the appended claims are intended to be interpreted as including both the preferred embodiments and all changes and modifications falling within the scope of the invention. Clearly, those skilled in the art can make various alterations and modifications to the invention without departing from its spirit and scope. Thus, if these modifications and modifications of the invention fall within the scope of the claims and their equivalents, the invention is also intended to include these modifications and modifications.

Claims

1. A Hamiltonian copy exchange dynamics method for self-assembled systems, characterized in that, include: S1, determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method; S2, based on the Hamiltonian replica exchange dynamics method, implement different Hamiltonian replica exchange dynamics transformations for the self-assembled system using either the REST method or the REST2 method; The Hamiltonian replica exchange dynamics transformation of the self-assembled system is implemented for the REST method based on the Hamiltonian replica exchange dynamics method, including: S21, determine the temperature of the m-th replica as T. m The biased potential function in this copy is: Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy; S22, determine the scaling factor λ in different replicas; and respectively for the interaction E pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the following equation (4) as the potential function under temperature T0: λ is greater than or equal to 0, and f(λ) is S23 transforms each replica in REST into a simulation at the same temperature T0, with the only difference being the biased potential function.

2. The Hamiltonian replica exchange dynamics method for self-assembled systems according to claim 1, characterized in that, The proportionality coefficient λ is a factor. Therefore, equation (4) is transformed into the following potential function: Therefore, the m-th replica can also be regarded as a system simulated by Equation (2) as the potential function under the condition of temperature T0.

3. The Hamiltonian replica exchange dynamics method for self-assembled systems according to claim 1, characterized in that, The Hamiltonian replica exchange dynamics transformation of the self-assembled system for the REST2 method is implemented based on the Hamiltonian replica exchange dynamics method, including: S21', determine the temperature of the m-th replica as T. m The biased potential function in this copy is: Where E pp E pw and E ww These represent solute-solute, solute-water, and water-water interactions, respectively. β equals 1 / kT, where k is the Boltzmann constant and T is the temperature; β0 equals 1 / kT0, where k is the Boltzmann constant and T0 is the initial temperature; β m equals 1 / kT m k is the Boltzmann constant, T m The temperature of the m-th copy; S22', determine the scaling factor λ in different copies; respectively for the interaction E pp and E pw Multiplying by the scaling factor λ and the scaling factor f(λ), the m-th replica is transformed into a system simulated by the potential function of equation (5) at temperature T0: S23' transforms each replica in REST2 into a simulation at the same temperature T0, with the only difference being the biased potential function.

4. The Hamiltonian replica exchange dynamics method for self-assembled systems according to claim 3, characterized in that, λ is greater than or equal to 0.

5. The Hamiltonian replica exchange dynamics method for self-assembled systems according to claim 4, characterized in that, The f(λ) is 6. The Hamiltonian replica exchange dynamics method for self-assembled systems according to claim 5, characterized in that, The proportionality coefficient λ is a factor. Therefore, equation (5) is transformed into the following potential function:

7. A Hamiltonian replica swapping dynamics system for self-assembly systems, used to implement the method according to any one of claims 1-6, characterized in that, include: The method determination module (101) is used to determine whether the Hamiltonian replica exchange dynamics method is the REST method or the REST2 method; The dynamic transformation module (102) is used to implement different Hamiltonian replica exchange dynamic transformations of the self-assembled system according to the Hamiltonian replica exchange dynamics method for the REST method or the REST2 method.

Images