Molecular dynamics-based coating surface hydration layer parallel analysis method and system
By employing a parallel molecular dynamics analysis method, the problem of low efficiency in hydration layer analysis during molecular dynamics simulations was solved. This method enables parallel computation of multiple attributes and result integration, thereby improving computational efficiency and the comprehensiveness of the analysis.
Patent Information
- Application Number
- CN202511994434.2
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-12-26
- Publication Date
- 2026-04-21
AI Technical Summary
Existing technologies are computationally inefficient in molecular dynamics simulations, struggle to handle hydration layer analysis of long-duration, multi-frame trajectories, and have fragmented analysis functions with a lack of a unified interface, making it impossible to achieve parallel computation of multiple attributes and integration of results.
A parallel analysis method for the hydration layer on the coating surface based on molecular dynamics is adopted. An analysis script is generated by trajectory segmentation, and multiple segmented trajectories are executed in parallel. By combining background and foreground processes, data such as the radial distribution function, orientation distribution function, and hydrogen bond structure of the hydration layer are obtained, and an analysis report is generated.
It enables comprehensive quantitative analysis of the structure, dynamics, and stability of the hydration layer, improves computational efficiency, avoids errors caused by switching between multiple tools, and achieves parallel computation of multiple attributes and result integration.
Smart Images

Figure CN121905310A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of molecular dynamics simulation technology, specifically relating to a parallel analysis method and system for hydration layers on coating surfaces based on molecular dynamics. Background Technology
[0002] The statements in this section are merely background information related to the present invention and do not necessarily constitute prior art.
[0003] Hydrophilic materials have attracted much attention due to their unique interfacial properties. At the molecular level, these materials form a stable hydration layer through surface hydrogen bonds and Coulomb interactions. This structural arrangement lays the foundation for their various superior properties in aqueous solutions, including inhibiting protein adsorption, enhancing interfacial lubrication, and improving protein stability. Currently, the presence of a hydration layer on the surface of polymer brush coatings can be determined experimentally. However, specific information about the hydration layer on the surface of polymer brush coatings, such as density, water molecule orientation, and hydrogen bond distribution, cannot be determined using microscopic observation methods such as nuclear magnetic resonance spectroscopy, X-ray diffraction, infrared spectroscopy, and atomic force microscopy.
[0004] Molecular dynamics simulations have become an important tool for studying the hydration behavior of material surfaces. By combining Newtonian mechanics with carefully parameterized interatomic potential energy, molecular dynamics provides atomic-scale insights into biomolecular systems, achieving sub-nanometer spatial resolution and femtosecond temporal resolution, surpassing the capabilities of most experimental techniques. Continuous optimization of force field parameters further ensures that simulation results closely match physical reality.
[0005] While some tools exist for calculating partial hydration layer parameters, they suffer from low computational efficiency and difficulty in handling long-duration, multi-frame trajectories. Traditional serial analysis methods require processing each frame of trajectory data sequentially. For large-scale, long-term simulated trajectories, this processing method is extremely time-consuming, becoming a major factor restricting research progress. Furthermore, the methods are fragmented and lack a unified interface; a comprehensive characterization of the hydration layer requires consideration from multiple dimensions (at least including structure, orientation, interaction, dynamics, and stability).
[0006] Currently, these analytical functions are scattered across different scripts, plugins, or software commands, requiring the manual running of multiple independent analytical processes. This is not only cumbersome, but more seriously, because these analyses may be based on different parameter settings, atom selections, or algorithm implementations, the final data cannot achieve multi-attribute parallel computation and result integration.
[0007] Therefore, there is an urgent need for a technical solution that can overcome the above-mentioned shortcomings, perform efficient parallel processing of large-scale molecular dynamics trajectories, integrate the calculation of multiple key properties of hydration layers, and provide a standardized and repeatable analysis process, thereby providing a powerful and reliable tool for interfacial hydration research in the fields of materials science and biophysics. Summary of the Invention
[0008] To address the aforementioned issues, this invention proposes a parallel analysis method and system for hydration layers on coating surfaces based on molecular dynamics. This method analyzes and performs parallel calculations on the structure of the hydration layer, enabling a comprehensive quantitative analysis of the structure, dynamics, and stability of the hydration layer. This solves the problems of low computational efficiency, fragmented functions, and insufficient automation in existing technologies.
[0009] According to some embodiments, the first aspect of the present invention provides a parallel analysis method for hydration layers on coating surfaces based on molecular dynamics, employing the following technical solution: A parallel analysis method for hydration layers on coating surfaces based on molecular dynamics includes: Obtain a molecular dynamics model of the surface hydration layer; The obtained molecular dynamics model is segmented into several segments. Based on the obtained segmentation trajectories, an analysis script for the segmentation trajectories is automatically generated; The analysis scripts for several segmented trajectories are executed in parallel. During the parallel execution, the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory are obtained. The spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer are analyzed to complete the parallel analysis of the coating surface hydration layer based on molecular dynamics.
[0010] As a further technical limitation, during the trajectory segmentation process, the obtained molecular dynamics model archive file is read to obtain the trajectory archive, and the obtained trajectory archive is segmented by the trajectory segmentation function to obtain several trajectory segmentation files, that is, several segmented trajectories.
[0011] As a further technical limitation, in the process of executing the analysis scripts of several segmented trajectories in parallel, a parallel computing method combining background processes and foreground processes is adopted, that is, a step-by-step execution strategy is adopted in conjunction with the task start prompt information.
[0012] Furthermore, the step-by-step execution strategy adopted is to set the first (n-1) computational tasks to run in the background, with each task executing its corresponding analysis script through an independent molecular dynamics visualization tool process running in the background, in order to avoid interfering with the main process.
[0013] As a further technical limitation, the average values of the obtained radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time and water binding site data are statistically averaged. An analysis report of the surface hydration layer is generated by combining the statistical average values. Based on the generated analysis report of the surface hydration layer, a parallel analysis of the coating surface hydration layer based on molecular dynamics is completed.
[0014] As a further technical constraint, before executing the analysis scripts for several segmented trajectories in parallel, the frame range, number of segments, step size, atomic selection statements, and radius of action in the obtained analysis scripts are validated for legality.
[0015] According to some embodiments, a second aspect of the present invention provides a parallel analysis system for hydration layers on coating surfaces based on molecular dynamics, employing the following technical solution: A parallel analysis system for hydration layers on coating surfaces based on molecular dynamics, comprising: The acquisition module is configured to acquire a molecular dynamics model of the surface hydration layer. The segmentation module is configured to segment the acquired molecular dynamics model into several segmented trajectories. The generation module is configured to automatically generate an analysis script for the segmented trajectories based on the obtained segmented trajectories. The analysis module is configured to execute analysis scripts for several segmented trajectories in parallel. During parallel execution, it acquires the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory. It analyzes the spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer, thus completing the parallel analysis of the coating surface hydration layer based on molecular dynamics.
[0016] According to some embodiments, a third aspect of the present invention provides a computer-readable storage medium, employing the following technical solution: A computer-readable storage medium having a program stored thereon, which, when executed by a processor, implements the steps in the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in the first aspect of the present invention.
[0017] According to some embodiments, the fourth aspect of the present invention provides an electronic device, which adopts the following technical solution: An electronic device includes a memory, a processor, and a program stored in the memory and running on the processor, wherein the processor executes the program to implement the steps in the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in the first aspect of the present invention.
[0018] According to some embodiments, the fifth aspect of the present invention provides a computer program product, which adopts the following technical solution: A computer program product includes software code, wherein the program in the software code performs the steps in the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in the first aspect of the present invention.
[0019] Compared with the prior art, the beneficial effects of the present invention are as follows: This invention integrates multiple core hydration layer property analyses into a unified framework, allowing users to initiate all calculations with a single command, avoiding errors caused by switching between multiple tools. It analyzes and performs parallel calculations on the surface hydration layer structure, achieving comprehensive quantitative analysis of the hydration layer structure, dynamics, and stability, thereby solving the problems of low computational efficiency, fragmented functions, and insufficient automation in existing technologies. Attached Figure Description
[0020] The accompanying drawings, which form part of this embodiment, are used to provide a further understanding of this embodiment. The illustrative embodiments and their descriptions are used to explain this embodiment and do not constitute an improper limitation of this embodiment.
[0021] Figure 1 This is a flowchart of the parallel analysis method for the hydration layer on the coating surface based on molecular dynamics in Embodiment 1 of the present invention; Figure 2 This is a schematic diagram of the molecular dynamics model of the monomer in Embodiment 1 of the present invention; Figure 3 These are schematic diagrams of the hydration peak structures of different systems in Embodiment 1 of the present invention; Figure 4 This is a schematic diagram of the orientation distribution function curves of different systems in Embodiment 1 of the present invention; Figure 5 This is a schematic diagram showing the hydrogen bond number distribution of different systems in Embodiment 1 of the present invention; Figure 6 This is a schematic diagram of the mean square displacement curves of different systems in Embodiment 1 of the present invention; Figure 7 This is a schematic diagram of the average residence time curves for different systems in Embodiment 1 of the present invention; Figure 8 This is a schematic diagram of the stable hydrogen bond curves of different systems in Embodiment 1 of the present invention; Figure 9 This is a structural block diagram of the parallel analysis system for the hydration layer of the coating surface based on molecular dynamics in Embodiment 2 of the present invention. Detailed Implementation
[0022] The present invention will be further described below with reference to the accompanying drawings and embodiments.
[0023] It should be noted that the following detailed descriptions are exemplary and intended to provide further illustration of the invention. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.
[0024] It should be noted that the terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the scope of exemplary embodiments according to the invention. As used herein, the singular form is intended to include the plural form as well, unless the context clearly indicates otherwise. Furthermore, it should be understood that when the terms "comprising" and / or "including" are used in this specification, they indicate the presence of features, steps, operations, devices, components, and / or combinations thereof.
[0025] In this invention, terms such as "upper," "lower," "left," "right," "front," "back," "vertical," "horizontal," "side," and "bottom" indicate the orientation or positional relationship based on the orientation or positional relationship shown in the accompanying drawings. These terms are used only to facilitate the description of the structural relationships of the various components or elements of this invention and do not specifically refer to any component or element in this invention. They should not be construed as limiting the invention.
[0026] In this invention, terms such as "fixed connection," "connected," and "linked" should be interpreted broadly, indicating a fixed connection, an integral connection, or a detachable connection; a direct connection or an indirect connection through an intermediate medium. Those skilled in the art can determine the specific meaning of these terms in this invention based on the specific circumstances, and they should not be construed as limitations on the invention.
[0027] Where there is no conflict, the embodiments and features in the embodiments of the present invention can be combined with each other.
[0028] Example 1 Embodiment 1 of this invention introduces a parallel analysis method for the hydration layer on the coating surface based on molecular dynamics.
[0029] like Figure 1 The method shown is a parallel analysis method for hydration layers on coating surfaces based on molecular dynamics, comprising: Obtain a molecular dynamics model of the surface hydration layer; The obtained molecular dynamics model is segmented into several segments. Based on the obtained segmentation trajectories, an analysis script for the segmentation trajectories is automatically generated; The analysis scripts for several segmented trajectories are executed in parallel. During the parallel execution, the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory are obtained. The spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer are analyzed to complete the parallel analysis of the coating surface hydration layer based on molecular dynamics.
[0030] Visual Molecular Dynamics (VMD) is a 3D visualization tool for modeling and analyzing biomolecular systems. It supports the display and animation of structures such as proteins and nucleic acids, and can be used to import PDB files, rotate from multiple perspectives (left mouse button), zoom (scroll wheel), and switch display modes (such as CPK ball-and-stick models). It is applied to molecular dynamics simulation trajectory analysis.
[0031] This embodiment utilizes a pseudo-parallel approach to rapidly calculate hydration layer structures in molecular dynamics simulations using VMD. By calculating and analyzing the surface radial distribution function (RDF), orientation distribution function (ODF), hydrogen bond structure (HB), mean square displacement (MSD), mean residence time (MRT), and water binding sites (WBS) of hydration layers formed by different polymer brushes, the true functionality and performance of this hydration layer calculator are evaluated. Since there are various polymer brush coating materials, but the molecular dynamics modeling and simulation process is consistent, this embodiment focuses on two such materials. Figure 2 All-atom molecular dynamics simulations were performed on common monomers shown in the figure—ethylene glycol (represented by PEG) and trimethylamine N-oxide (represented by TMAO).
[0032] This embodiment employs molecular dynamics simulation to study the polymer brush coating system. Force field parameters are generated using the CGenFF tool, compatible with the CHARMM36 force field. The water model is TIP3P, and CMAP correction is introduced to optimize the dihedral angle description of the main chain. The model is constructed using VMD, and the simulation is performed in NAMD 3.0 software. The system is set to periodic boundary conditions, with a simulation box dimension of 100 Å (length, width, and height), a polymer brush placed at the center, and water molecules added as solvent. A pure water system is also established as a control. To maintain system stability, positional constraints are applied to the terminal carbon atoms of the monomers and the oxygen atoms in the water molecules.
[0033] The simulation first performed energy minimization, followed by a 10 ns NPT equilibration at 1 atm, and then a 40 ns production simulation under NVT conditions. The time step was 2 fs. Temperature control was achieved using the Langevin method (298.15 K), and pressure control was achieved using the Langevin piston. Long-range electrostatic interactions were handled using the PME method, and van der Waals interactions were double-trunculated (10 / 12 Å). Trajectories were saved every 2 ps for analysis.
[0034] In this embodiment, the user needs to open the VMD software and enter the Tk Console to activate the Tcl / Tk script interface. It is essential to ensure that the Tcl script file, molecular structure file, and trajectory file (specifically in .psf and .dcd formats) for analysis are saved in the current working directory of VMD to ensure proper file access permissions. To start the analysis, the user simply needs to enter the source file name of the Tcl script in the Tk Console to load the script; the tool can be executed using the following command format: evaluate psffile dcdfile -ref $value option1 $value1 option2 $value2... Users need to input the values of required parameters such as psffile, dcdfile, and -ref. All other parameters will be automatically assigned values or can be manually adjusted according to actual needs. The usage and requirements of all available parameters are as follows: psffile, the structure file format should be .psf; dcdfile, the coordinate file format should be .dcd; -first, specifies the starting frame, which should be an integer, and the default is 0, i.e., the first frame. -last; specifies the end frame, should be an integer, defaults to -1 (last frame); -n, the number of segments in the .dcd file, defaults to 8; -s, the step size between calculations, defaults to 1; -ref, the reference atom for forming the hydration layer; -cal, the atom being calculated, usually H2. O atoms in O; default value, bond type is OT; -r, radius of the hydration layer, default value is 3.5; -RDF, option to calculate the atomic radial distribution function, 1 indicates calculation, 0 indicates skip, default is 0; -ODF, option to calculate the atomic orientation distribution function, 1 indicates calculation, 0 indicates skip, default is 0; -HB, option to calculate the interatomic hydrogen bond structure, 1 indicates calculation, 0 indicates skip, default is 0; -MSD, option to calculate the atomic mean square displacement, 1 indicates calculation, 0 indicates skip, default is 0; -MRT, option to calculate the atomic mean residence time, 1 indicates calculation, 0 indicates skip, default is 0; -WBS, option to calculate the atomic water-binding sites, 1 indicates calculation, 0 indicates skip, default is 0.
[0035] The simulation process of each calculation module is mainly divided into two stages. In the first stage, the corresponding module's calculation function is called based on the input parameters or default parameters. In the second stage, the segmented calculation results are integrated and the results are output to Data_. $value .dat file.
[0036] Module 1: Atomic Radial Distribution Function (RDF) Calculation, specifically: RDF analysis begins with the measure_RDF procedure, which includes receiving the output filename, reference atom selection statement, computation atom selection statement, and optional parameters such as the spherical shell distance resolution (dr, default 0.1 Å) and the maximum distance from the reference atom (rmax, default 10.0 Å). The procedure first initializes the computational environment, including opening the output file, obtaining the total number of simulation frames, creating the selection sets of reference and computation atoms, and counting the number of reference and computation atoms.
[0037] The computational domain is divided into multiple distance intervals (i.e., spherical shells of different sizes) (bins), the number of which is determined by rmax and dr, i.e., the number of spherical shells variable numbin = int(rmax / dr). Each interval corresponds to an intermediate distance value variable rmid =(j + 0.5) * dr, where j is the interval index. Simultaneously, the system calculates the volume V of the simulation box and the average density avgp_cal of atoms for subsequent integration processing.
[0038] While traversing each frame of the simulated trajectory, the process calculates the distribution of computed atoms around each reference atom. For each distance interval, the VMD's `atomselect` command, combined with the `pbwithin` condition, selects computed atoms within the range of the variable `rtmp` (the lower bound of the current distance interval) and the variable `rtmp+dr`, respectively. By calculating the difference in the number of atoms in the two selection sets, the number of newly added atoms in that interval is obtained and accumulated into the corresponding `rmid` count array `num`. This step ensures the accuracy of the atom count within each distance interval while considering periodic boundary conditions. After traversing all frames and all reference atoms, the process performs post-processing on the statistical data. For each distance interval, the actual atom count variable `dn` (i.e., `num(rmid)` divided by the total number of frames and the number of reference atoms) is calculated, and the corresponding spherical shell volume variable `dV` is calculated. Finally, the RDF value is calculated as `dn / dV` divided by the average density of computed atoms `avgp_cal`, thus obtaining the radial distribution function `g(r)`. The results are written to the output file in the form of distance `rmid` and the corresponding `g(r)` value.
[0039] After completing the RDF calculation for each sub-trajectory, the `integrate_RDF` procedure integrates the results of multiple parallel tasks and performs data processing. This procedure receives the number of subtasks, `n`, as input and initializes the storage data structure. For each subtask's result file (e.g., `Data_RDF_1.dat`, `Data_RDF_2.dat`), the procedure reads the distance and g(r) values from each row. On the first read, the distance values are stored as the baseline column; subsequent subtask data have their g(r) values extracted row-by-row and stored in a list. During the integration phase, for each distance point, the average of the corresponding g(r) values from all subtasks is calculated to reduce potential random errors and improve statistical reliability. Finally, the averaged distance and g(r) values are written to the main output file of the radial distribution function, `Data_RDF.dat`, forming the final RDF result file.
[0040] A distinct hydration peak was observed in all three systems in the example (bulk water, PEG, and TMAO). The main peak positions were approximately 3 Å, 3 Å, and 2.8 / 5.8 / 6.2 Å, respectively (e.g., ...). Figure 3 This confirms the existence of at least one well-defined hydration layer structure. TMAO shows three peaks, but the first peak is clearly more prominent, indicating a highly ordered first shell and a more diffuse distribution of the second and third shells. The relative intensities of the hydration peaks are ordered as water > TMAO > PEG, consistent with its known hydrophilic characteristics. Due to the weaker intensity of the second peak, the hydration layer boundary is set at 3.5 Å in this embodiment.
[0041] Module 2: Atomic Orientation Distribution Function (ODF) Calculation, specifically: The orientation distribution function (ODF) is an important property characterizing the orientation of water molecules relative to a reference atom, used to analyze the distribution of the water molecule dipole moment around the reference atom. ODF analysis begins with the `measure_ODF` procedure, which receives the output filename, reference atom selection statement, calculation atom selection statement, and the radius of influence (rmax) parameter. The procedure first initializes the computational environment, including opening the output file, obtaining the total number of simulation frames, and creating the selection set of reference atoms. The system initializes a frequency array `frequency` to record the occurrence count of each integer angle value from 1 to 180 degrees, and simultaneously initializes the total count variable `totalnum` to zero.
[0042] When traversing each frame of the simulated trajectory, the process first updates the coordinate information of the reference atom. For each reference atom, the VMD's `atomselect` command, combined with the `pbwithin` condition, selects water molecules (computation atoms) within the `rmax` range surrounding the reference atom. For each selected water molecule, the process further identifies its oxygen atom (OT) and two hydrogen atoms (HT), obtaining the hydrogen atom index through the bonding information of the oxygen atom. The system creates a corresponding atom selection set, including the oxygen atom, the two hydrogen atoms, and the reference atom, and updates their coordinate information. The orientation vector of the water molecule is calculated, i.e., the vector from the oxygen atom (the positive dipole of the water molecule) to the center of the two hydrogen atoms (the negative dipole of the water molecule), and the vector from the reference atom to the oxygen atom. The cosine of the angle between these two vectors is calculated using the vector dot product formula, and then the radian value is converted to an angle value using the inverse cosine function. Finally, the integer angle value variable `Ang` is obtained using the round-up function. The frequency count variable `frequency(Ang)` of this angle value is incremented by 1, and the total count variable `totalnum` is also incremented by 1.
[0043] After traversing all frames and all reference atoms, the process writes the total count (totalnum) to the first line of the output file, followed by the frequency values corresponding to each angle from 1 to 180 degrees in sequence. This output format facilitates subsequent data processing and integration.
[0044] After calculating the ODF for each sub-trajectory, the `integrate_ODF` procedure integrates the results of multiple parallel tasks and processes the data. This procedure receives the number of sub-task result files, `n`, as input and initializes the data storage array `data`. For each sub-task result file, the procedure reads the total count from the first line and the angular frequency values from the subsequent 180 lines. The total count from the first line is added to the overall total count, while the frequency values for each angle are added to the data array according to their corresponding positions. During the integration phase, for each angle (from 1° to 180°), the sum of the corresponding frequency values from all sub-task result files is calculated, then divided by the overall total count and multiplied by 100 to obtain the percentage of occurrence for that angle. Finally, the overall total count is written to the first line of the output file `Data_ODF.dat`, followed by the percentage values for each angle from 1° to 180°.
[0045] In all three systems in this embodiment, the highest ODF peak appears around 130° (e.g., Figure 4The presence of an ordered hydration layer is indicated. Further analysis revealed that TMAO exhibits a particularly sharp peak, reflecting the highly ordered arrangement of water molecules around the monomer, suggesting excellent hydrophilicity. In contrast, the peaks of other materials are significantly broader and less intense. Notably, both bulk water and PEG show a secondary peak near 50°, which can be interpreted as evidence of two different hydrogen bonding modes. The difference in peak intensity can be explained by hydrogen bond configurations: when water molecules act as hydrogen bond donors and bind to the reference atom, the rotational degrees of freedom of the relevant hydrogen atoms are restricted, resulting in a more pronounced ODF peak. Conversely, when water molecules act as hydrogen bond acceptors, their hydrogen atoms retain greater rotational degrees of freedom, leading to a broader peak and lower intensity.
[0046] Module 3 Hydrogen Bond Structure Analysis: Specifically, hydrogen bond analysis begins with the measure_HB procedure (this tool defines hydrogen bond structure by default as: donor). The receptor distance is less than 3.5 angstroms and the receptor The donor-hydrogen bond angle is less than 30 degrees. This process receives the output filename, reference atom selection statement, calculation atom selection statement, interaction radius (rmax), and angle threshold (amax, default 30°). The process first initializes the computational environment, including opening the output file, obtaining the total number of simulation frames, and creating a selection set of reference atoms. The process converts the angle to radians and calculates the cosine value as the geometric criterion for hydrogen bond formation. Additionally, the system writes the total number of frames and the number of reference atoms to the first line of the output file.
[0047] When traversing each frame of the simulated trajectory, the process analyzes hydrogen bond formation in two ways. First, when the reference atom acts as a hydrogen bond acceptor, water molecules within a radius of rmax around each reference atom are selected. For each water molecule, its oxygen atom and two hydrogen atoms are obtained, and the vector from the oxygen atom to each hydrogen atom and the vector from the reference atom to the oxygen atom are calculated. The cosine of the angle between the relevant vectors of each hydrogen atom is calculated using the vector dot product. If any cosine value is greater than a preset threshold, a hydrogen bond is determined to have formed, and the hydrogen bond type is marked as "1" (acceptor type). The hydrogen bond distance and angle values are recorded to the output file. Second, when the reference atom acts as a hydrogen bond donor, its bonded atoms are checked for each reference atom, and hydrogen atoms with a mass rounded to 1 amu are selected. For each bonded hydrogen atom, water molecules and oxygen atoms within a radius of rmax around the reference atom are selected, and the vector from the reference atom to the bonded hydrogen atom and the vector from the reference atom to the water-oxygen atom are calculated. The cosine of the included angle is calculated by vector dot product. If it is greater than the threshold, a hydrogen bond is determined to be formed, and the hydrogen bond type is marked as "0" (donor type). At the same time, the corresponding distance and angle values are recorded.
[0048] After completing the hydrogen bond analysis for each sub-trajectory, the `integrate_HB` process is responsible for integrating the results of multiple parallel tasks and performing data processing. This process receives the number of result files `n` from the sub-tasks, as well as optional parameters such as distance range, angle range, and analysis precision. The process initializes the total number of frames `nframe` and the number of reference atoms `numref`; initializes the counting arrays `count(0)` and `count(1)`, which are used to count the total number of donor and acceptor hydrogen bonds, respectively; calculates the number of distance intervals `num_dist_bins` and the number of angle intervals `num_angle_bins`; and initializes the distribution arrays `dist_bins(0),` `dist_bins(1)`, `angle_bins(0),` and `angle_bins(1)`, which are used to record the distribution of the two types of hydrogen bonds in the distance and angle intervals. For each sub-task's result file, the process first reads the number of frames and reference atoms from the first line and adds them to the total count. Then, it reads the hydrogen bond records line by line, counts the number of acceptor and donor hydrogen bonds according to the type identifier, and simultaneously counts their distribution in the distance and angle intervals. In the integration phase, the process first calculates the normalization factor, which is the product of the total number of frames and the number of reference atoms. Based on this, the probabilities of acceptor-type hydrogen bonds and donor-type hydrogen bonds are calculated (acceptor-type hydrogen bond probability = count(1) / normalization factor; donor-type hydrogen bond probability = count(0) / normalization factor). These two probabilities correspond to the average number of hydrogen bonds formed when each reference atom acts as an acceptor and donor, respectively. For distance and angular distributions, the percentage of hydrogen bonds in each interval relative to the total number of hydrogen bonds of that type is calculated, thus obtaining their distribution curves. The final output file contains hydrogen bond probability data as well as detailed distance and angular distribution information.
[0049] In all three systems of this embodiment, the number of hydrogen bonds formed by the reference atom is as follows: Figure 5 As shown in the figure. Analysis indicates that the ratio of acceptor to donor interactions in bulk water is essentially balanced, which is entirely consistent with expectations. TMAO, as an acceptor, forms an average of 1.9 hydrogen bonds, while the donor value is 0. This directly corresponds to its molecular structure: TMAO lacks hydrogen atoms covalently bonded to oxygen or nitrogen (its three hydrogen atoms are all attached to carbon), thus it cannot act as a donor for traditional hydrogen bonds. This characteristic makes TMAO primarily a strong hydrogen bond acceptor in solution. When it binds to water molecules, it fixes the donor hydrogen atoms in the water molecules, thereby interfering with the original hydrogen bond connections between these water molecules and altering the dynamics and orderliness of the local hydration structure. PEG, as an acceptor, forms an average of about 1.2 hydrogen bonds, while as a donor, it forms an average of 0.79 hydrogen bonds. The acceptor-donor interaction is slightly unbalanced, and both are less than those in bulk water, indicating that the hydrogen bond interaction between PEG and water is weaker in both strength and quantity than the network formed by water itself. This allows it to dissolve in water without causing drastic disturbance to the water phase structure.
[0050] Module 4, Atom Mean Square Displacement Analysis, specifically begins with the `measure_MSD` procedure. This procedure receives the output filename, reference atom selection statement, calculation atom selection statement, and radius of action parameter. The procedure first initializes the output file and obtains the total number of simulation frames, creating a selection set of reference atoms. For each frame from the initial time to 100 frames (total frames minus 100), the system selects the reference atom for the current frame and analyzes each water molecule within its radius of action. For each water molecule, the procedure records the squared change in distance relative to the reference atom over the time span from the current frame to the next 50 frames. In the implementation, the system first obtains the initial distance between the water molecule and the reference atom as a baseline value. Then, for each subsequent frame, it calculates the squared difference between the current distance and the initial distance and organizes these squared differences into a data sequence in chronological order. The entire process iterates through all frames and all reference atoms that meet the criteria.
[0051] After calculating the mean square displacement of each sub-trajectory and outputting the result files for each subtask, the results of all parallel subtasks are integrated using the `integrate_MSD` procedure. This procedure receives the number of subtasks as input and initializes the data storage array. During the data reading phase, the procedure iterates through the result files of all subtasks. For each mean square displacement record in each file, the system stores these values column-wise into different arrays. All square displacement values at each time point are collected into the corresponding data list. The integrated result first outputs the total number of valid water molecule trajectories, and then calculates the arithmetic mean of all square displacement values at each time point to eliminate random fluctuations of individual water molecules.
[0052] In all three systems of this embodiment, the resulting MSD curves are as follows: Figure 6As shown in the figure, the curve representing TMAO exhibits the highest growth rate and the largest displacement amplitude, indicating that under unconstrained interface conditions, water molecules demonstrate rapid free diffusion within a simulated 50 ps, with a mobility even exceeding that of bulk water. This is attributed to TMAO's unique amphiphilic structure: its strongly polar NO group has a dipole interaction with water, but its hydrophobic shell composed of three methyl groups prevents the formation of persistent, network-like hydrogen bonds with water. This characteristic makes it less susceptible to being tightly bound by the hydration layer, thus exhibiting an exceptionally high mobility. This phenomenon is consistent with TMAO's known function in biology, enabling rapid movement inside and outside cells. The PEG curve in the middle presents a stark contrast, with a much lower MSD value than TMAO, indicating extremely slow growth. This directly confirms the highly restricted movement of PEG segments (or its monomers). Strong hydrogen bonding, conformational entropy constraints of the polymer chain, and the potential formation of entangled networks collectively "anchor" it near the interface, resulting in extremely low diffusion capacity. Its movement pattern is mainly local oscillation and relaxation, rather than long-range translation.
[0053] Module 5, Atom Average Residence Time Analysis, specifically begins with the `measure_MRT` procedure. This procedure receives the output filename, reference atom selection statement, calculation atom selection statement, and radius of influence parameter. The procedure first initializes the output file and obtains the total number of simulation frames, creating a selection set of reference atoms. For each frame from the first frame to the total frame minus one hundred, the system selects the reference atom for the current frame and dynamically tracks the behavior of water molecules within its radius of influence. For each water molecule, the procedure establishes a monitoring window of fifty-one time points, verifying the initial state of the water molecule starting from the frame preceding the current frame. The system first checks whether the water molecule was already within the radius of influence of the reference atom in the previous frame; if so, it skips the molecule to avoid duplicate counting. For newly entering water molecules within the radius of influence, the procedure initializes a list of length 51, with all elements initialized to 1; it also sets a flag `found` to 0, indicating the hypothetical state that the water molecule has not left and continues to remain within the radius of influence. Subsequently, the system detects the change in distance between water molecules and reference atoms frame by frame from the current frame to fifty frames later. When a water molecule first leaves the effective radius, the system updates the marker value to zero for that moment and all subsequent time points, and sets `found` to 1 to record the water molecule's departure. The entire process iterates through all frames that meet the criteria and all reference atoms to ensure sufficient statistical sampling.
[0054] After calculating the average dwell time of each sub-trajectory and outputting the result file of the subtask, the `integrate_MRT` procedure is responsible for integrating the result files of all parallel subtasks. This procedure receives the number of subtasks `n` as input and initializes the column data summation array `column_sums` and the row counter variable `total_lines`. During the data reading phase, the procedure iterates through the result files of all subtasks. For each marker sequence record in each file, the system processes the data column by column. For each time point column, the system accumulates all values in that column and simultaneously counts the total number of valid trajectories. This processing method corresponds to different water molecule trajectories vertically and fifty-one consecutive frames horizontally. In the integration phase, the total number of trajectories involved in the statistics is first output, and then the arithmetic mean of all marker values at that time point is calculated for the corresponding column of each frame. This averaging process converts the binary markers into a continuous probability curve, and the average value at each time point represents the probability that a water molecule remains within its effective range from the initial time to that time point. The final output file contains the total number of trajectories and the average retention probability corresponding to the fifty-one time points, forming a complete time curve of water molecule retention probability.
[0055] In all three systems of this embodiment, such as Figure 7 As shown, the shorter MRT indicates a highly efficient water exchange rate between TMAO and PEG and the bulk phase. This dynamic conversion mechanism is crucial for maintaining effective lubrication, which makes charged polymer materials unique in antifouling surfaces, perfectly consistent with their expected relative hydrophilicity.
[0056] Module 6, Atom-Binding Water Site Analysis, specifically begins with the `measure_WBS` procedure. This procedure receives the output filename, reference atom selection statement, calculation atom selection statement, interaction radius parameter, and angle threshold (default 30°). The procedure first initializes the output file and obtains the total number of simulation frames, creates a selection set of reference atoms, and removes hydrogen atoms based on atomic mass. The analysis employs a hydrogen bond event triggering mechanism, traversing each non-hydrogen reference atom and each frame of simulation data. For each reference atom, the system selects water molecules within its surrounding interaction radius and analyzes hydrogen bond formation from both perspectives. When a reference atom acts as a hydrogen bond donor, the system uses the VMD built-in `measure_hbonds` command to detect whether a hydrogen bond has formed between the reference atom and the oxygen atom of a water molecule. This command uses both distance and angle for judgment. After detecting the presence of a hydrogen bond, the system further verifies whether the hydrogen bond already existed in the previous frame. If it does not exist, it is considered a new hydrogen bond event, and stability tracking is initiated. The system initializes a 51-bit list of markers, with all elements set to the number 1, representing the hypothetical state of continued hydrogen bond existence. Then, the hydrogen bond state is detected frame by frame within a range of 50 frames. When a hydrogen bond breaks for the first time, the marker value of that moment and all subsequent frames is updated to zero.
[0057] The same analytical workflow applies to scenarios where the reference atom acts as a hydrogen bond acceptor. Each valid hydrogen bond event is recorded with its reference atom index, water molecule index, starting frame, and complete label sequence. After calculating the water binding stability of each sub-trajectory, the `integrate_WBS` process integrates the results of all parallel tasks. This process receives the number of result files from the sub-tasks and the stability threshold parameter, and initializes the data structures for the average value mapping dictionary and the count mapping dictionary. During the data reading phase, the process iterates through all the result files from the sub-tasks, parsing each hydrogen bond event record. The system first checks if the first `tcut` time points of the label sequence are all 1s; hydrogen bond events meeting this condition are considered stable binding events, and their lifetime information is recorded in a dedicated result file. Simultaneously, the system accumulates the label columns of all hydrogen bond events for each reference atom according to the reference atom index data and counts the number of events, forming a stability result based on each reference atom. In the integration phase, the system calculates the average label value for each reference atom at each time point, and obtains the hydrogen bond presence probability by dividing the accumulated value by the number of events. The final output file contains the identifier of each reference atom and the corresponding average hydrogen bond presence probability value within a 51-frame range.
[0058] In all three systems of this embodiment, such as Figure 8As shown, during the simulation, bound water molecules on the solute surface formed stable hydrogen bonds with the O2 of PEG and the O25 of TMAO, and were therefore classified as bound water. This is mainly attributed to the high local charge density around the oxygen atoms. Under the influence of Coulomb forces, these atoms are more likely to form stable hydrogen bonds with water molecules. Therefore, these regions are more likely to become sites for bound water.
[0059] To evaluate the computational efficiency of the code, this embodiment uses two methods to perform ODF calculations on the same system: one is a script that does not split files, and the other is the HydroLayer calculator. The original single-threaded method only supports single-threaded operation, while the HydroLayer calculator can enable multi-threaded operation based on the specified parameter -n value, fully utilizing CPU efficiency and thus significantly improving computational efficiency. For example, when the parameter -n is set to 4, the HydroLayer calculator's computation time is only about one-quarter of that of the single-threaded script, fully demonstrating its significant advantage in computational performance.
[0060] This embodiment systematically solves the three major stability problems caused by parallelization itself—process out of control, resource contention, and result integrity—through a combination of central scheduling, controlled parallelism, and proactive output verification. By using trajectory segmentation and multi-process parallelism, the multi-core performance of the processor can be fully utilized, reducing computation time by several to tens of times, making ultra-long trajectory analysis feasible. The analysis of six core hydration layer attributes is integrated into a unified framework, allowing users to start all calculations with a single command, avoiding errors caused by switching between multiple tools. Through parameter verification and error handling mechanisms, process management and result monitoring processes are introduced to prevent overall analysis interruption due to the failure of a single task, improving system stability. Automatic result collection, verification, averaging, and output ensure data consistency and reduce the user's burden.
[0061] Example 2 Embodiment 2 of the present invention introduces a parallel analysis system for hydration layers on coating surfaces based on molecular dynamics.
[0062] like Figure 9 The system shown is a parallel analysis system for hydration layers on coating surfaces based on molecular dynamics, comprising: The acquisition module is configured to acquire a molecular dynamics model of the surface hydration layer. The segmentation module is configured to segment the acquired molecular dynamics model into several segmented trajectories. The generation module is configured to automatically generate an analysis script for the segmented trajectories based on the obtained segmented trajectories. The analysis module is configured to execute analysis scripts for several segmented trajectories in parallel. During parallel execution, it acquires the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory. It analyzes the spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer, thus completing the parallel analysis of the coating surface hydration layer based on molecular dynamics.
[0063] The detailed steps are the same as those of the parallel analysis method for the hydration layer on the coating surface based on molecular dynamics provided in Example 1, and will not be repeated here.
[0064] Example 3 Embodiment 3 of the present invention provides a computer-readable storage medium.
[0065] A computer-readable storage medium having a program stored thereon, which, when executed by a processor, implements the steps in the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in Embodiment 1 of the present invention.
[0066] The detailed steps are the same as those of the parallel analysis method for the hydration layer on the coating surface based on molecular dynamics provided in Example 1, and will not be repeated here.
[0067] Example 4 Embodiment 4 of the present invention provides an electronic device.
[0068] An electronic device includes a memory, a processor, and a program stored in the memory and running on the processor. When the processor executes the program, it implements the steps in the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in Embodiment 1 of the present invention.
[0069] The detailed steps are the same as those of the parallel analysis method for the hydration layer on the coating surface based on molecular dynamics provided in Example 1, and will not be repeated here.
[0070] Example 5 Embodiment 5 of the present invention provides a computer program product.
[0071] A computer program product includes software code, wherein the program in the software code performs the steps of the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in Embodiment 1 of the present invention.
[0072] The detailed steps are the same as those of the parallel analysis method for the hydration layer on the coating surface based on molecular dynamics provided in Example 1, and will not be repeated here.
[0073] Those skilled in the art will understand that embodiments of the present invention can be provided as methods, systems, or computer program products. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the present invention can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code. The solutions in the embodiments of the present invention can be implemented using various computer languages, such as the object-oriented programming language Java and the interpreted scripting language JavaScript.
[0074] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.
[0075] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.
[0076] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.
[0077] 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.
[0078] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from its spirit and scope. Therefore, if these modifications and variations fall within the scope of the claims of this invention and their equivalents, this invention also intends to include these modifications and variations.
[0079] The above description is merely a preferred embodiment of this practice and is not intended to limit the scope of this practice. Various modifications and variations can be made to this practice by those skilled in the art. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of this practice should be included within the protection scope of this practice.
Claims
1. A parallel analysis method for hydration layers on coating surfaces based on molecular dynamics, characterized in that, include: Obtain a molecular dynamics model of the surface hydration layer; The obtained molecular dynamics model is segmented into several segments. Based on the obtained segmentation trajectories, an analysis script for the segmentation trajectories is automatically generated; The analysis scripts for several segmented trajectories are executed in parallel. During the parallel execution, the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory are obtained. The spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer are analyzed to complete the parallel analysis of the coating surface hydration layer based on molecular dynamics.
2. The parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in claim 1, characterized in that, During the trajectory segmentation process, the obtained molecular dynamics model archive file is read to obtain the trajectory archive. The obtained trajectory archive is then segmented using the trajectory segmentation function to obtain several trajectory segmentation files, which are in turn several segmented trajectories.
3. The parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in claim 1, characterized in that, During the parallel execution of analysis scripts for several segmented trajectories, a parallel computing method combining background and foreground processes is adopted, that is, a step-by-step execution strategy is adopted in conjunction with task startup prompt information.
4. The parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in claim 3, characterized in that, The step-by-step execution strategy adopted is to set the first (n-1) calculation tasks to run in the background, and each task executes the corresponding analysis script through an independent molecular dynamics visualization tool process in the background, so as to avoid interfering with the main process.
5. The parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in claim 1, characterized in that, The average values of the obtained radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data are statistically analyzed. An analysis report of the surface hydration layer is generated based on the statistical average values. Parallel analysis of the coating surface hydration layer based on molecular dynamics is then performed based on the generated analysis report.
6. The parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in claim 1, characterized in that, Before executing the analysis scripts for several segmented trajectories in parallel, the frame range, number of segments, step size, atomic selection statements, and radius of action in the obtained analysis scripts are validated for legality.
7. A parallel analysis system for hydration layers on coating surfaces based on molecular dynamics, characterized in that, include: The acquisition module is configured to acquire a molecular dynamics model of the surface hydration layer. The segmentation module is configured to segment the acquired molecular dynamics model into several segmented trajectories. The generation module is configured to automatically generate an analysis script for the segmented trajectories based on the obtained segmented trajectories. The analysis module is configured to execute analysis scripts for several segmented trajectories in parallel. During parallel execution, it acquires the radial distribution function, orientation distribution function, hydrogen bond structure, mean square displacement, average residence time, and water binding site data of the surface hydration layer in each segmented trajectory. It analyzes the spatial structure of the surface hydration layer, the directionality of the interaction between water molecules and the surface, the strength of the hydrogen bond network, and the dynamic properties of the hydration layer, thus completing the parallel analysis of the coating surface hydration layer based on molecular dynamics.
8. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the program is executed by the processor, it implements the steps of the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in any one of claims 1-6.
9. An electronic device comprising a memory, a processor, and a computer program stored in the memory and running on the processor, characterized in that, When the processor executes the program, it implements the steps of the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in any one of claims 1-6.
10. A computer program product, comprising software code, characterized in that, The program in the software code executes the steps of the parallel analysis method for hydration layers on coating surfaces based on molecular dynamics as described in any one of claims 1-6.
