High-throughput design method for crystal structures based on molecular or ionic configurations

Abstract

The application provides a crystal structure high-throughput design method based on molecular or ion configuration, which comprises the following steps: obtaining candidate ion structure data and performing quantum chemistry calculation and electrostatic potential analysis, combining ion crystal density empirical formula calculation to predict the crystal density, and screening ion combinations with high predicted density; further generating multiple initial crystal structure models based on space group theory and Wyckoff position symmetry analysis, and realizing efficient optimization of the crystal structure through two-level optimization steps of pre-optimization by a machine learning potential model and fine optimization by first principles; and performing bonding state analysis, protonation phenomenon detection and energy and density data extraction on the optimized structure, and screening the crystal structure meeting the performance index, so as to solve the problem that it is difficult to efficiently screen and generate high-density metastable crystal structures.

Description

Technical Field

[0001] This application relates to the field of energetic materials technology, and in particular to a high-throughput design method for crystal structures based on molecular or ionic configurations. Background Technology

[0002] In the research of energetic materials and other molecular / ionic crystal materials, the microscopic crystal structure is usually composed of molecules or ionic groups, and its configuration space expands exponentially with the increase of the types and stacking modes of molecules / ions. Therefore, how to efficiently generate, search for, and optimize molecular / ionic crystal structures, especially how to selectively screen and design paired molecules / ions and their corresponding crystal structures to improve crystal packing density, has become an urgent problem to be solved.

[0003] For molecular crystal structure design, existing methods employ symmetry-guided enumeration strategies. These methods systematically coarse-grained the configuration space using Wyckoff bit encoding, combining grid generation with a two-stage screening process using machine learning / density functional theory to efficiently explore potential energy surfaces. These methods aim to search for steady-state crystal structures based on energy minimization.

[0004] However, the aforementioned methods are primarily geared towards neutral molecular systems and are difficult to apply to the design of charged ionic crystal structures. Furthermore, since the target structures are often in metastable states and frequently require more compact crystal packing and higher density, traditional energy-minimizing steady-state search strategies are ineffective in capturing such high-density metastable structures. Therefore, existing technologies suffer from insufficient applicability to charged ionic systems and difficulties in efficiently screening and generating high-density metastable crystal structures. Summary of the Invention

[0005] This application provides a high-throughput design method for crystal structures based on molecular or ionic configurations to solve the problem of difficulty in efficiently screening and generating high-density metastable crystal structures.

[0006] This application provides a high-throughput design method for crystal structures based on molecular or ionic configurations, the method comprising: Acquire candidate ion structure data, which includes structural information of cations and anions; Quantum chemical calculations are performed on the ions in the candidate ion structure data to obtain the ion geometry and wavefunction; Based on the wave function, an electrostatic potential analysis is performed on the ion to obtain the molecular surface electrostatic potential distribution parameters of the ion; Based on the molecular surface electrostatic potential distribution parameters and the preset empirical formula for ionic crystal density, the predicted crystal density of ionic crystals formed by different combinations of cations and anions is calculated. Based on a preset density threshold, target ion combinations with a predicted crystal density greater than the density threshold are selected from the different combinations of cations and anions. Based on space group theory and Wyckoff position symmetry analysis, the molecular symmetry axes of the cations and anions in the target ion combination are matched with the symmetry axes of the selected Wyckoff positions in the target space group, and all symmetry operations of the target space group are applied to generate multiple initial crystal structure models. The multiple initial crystal structure models were pre-optimized using a machine learning potential model to obtain the first optimized crystal structure model. A first structure evaluation is performed on the first optimized crystal structure model to identify the bonding states between ions and obtain the energy and density data of the first crystal. Based on the results of the first structure evaluation, a first candidate crystal structure model that meets the first preset energy and density index is selected from the first optimized crystal structure model, and the first candidate crystal structure model is finely optimized using the first principle calculation method to obtain the second optimized crystal structure model. A second structure evaluation is performed on the second optimized crystal structure model to identify interionic bonding states, detect protonation, and obtain crystal energy and density data. Based on the results of the second structure evaluation, crystal structures that meet the preset energy and density performance indicators are selected from the second optimized crystal structure model.

[0007] By acquiring candidate ion structure data and performing quantum chemical calculations and electrostatic potential analysis, and combining this with empirical formulas for ion crystal density to predict crystal density, high-predicted density ion combinations are screened out. Furthermore, multiple initial crystal structure models are generated based on space group theory and Wyckoff position symmetry analysis. A two-stage optimization process—pre-optimization using machine learning potential models and fine-tuning using first-principles calculations—achieves efficient optimization of the crystal structure. The optimized structures are then subjected to bonding state analysis, protonation detection, and energy and density data extraction to screen out crystal structures that meet performance indicators. This method helps alleviate the efficiency problem of screening and generating high-density metastable crystal structures and improves the design process for energetic materials and other multi-component functional materials.

[0008] Optionally, the step of performing quantum chemical calculations on the ions in the candidate ion structure data includes: Within the framework of density functional theory, hybrid functionals and double zeta-polarized basis sets are used to perform geometric optimization calculations and frequency calculations for all ions, and the corresponding wavefunction files are output.

[0009] When performing quantum chemical calculations, by using hybrid functionals and double ζ-valence layer polarized basis sets within the framework of density functional theory to optimize the geometry and calculate the frequency of all ions and output wavefunction files, a balance can be achieved between computational accuracy and resource consumption. This provides a reliable data foundation for subsequent wavefunction-based electrostatic potential analysis and density prediction, thereby improving the accuracy of ion screening and the robustness of subsequent design processes.

[0010] Optionally, the electrostatic potential analysis includes: Based on the wavefunction file, the characteristic parameters of the molecular surface electrostatic potential distribution of each ion are calculated and extracted using wavefunction analysis software. The characteristic parameters include the surface area and average value of the part of the ion with positive electrostatic potential, and the surface area and average value of the part of the ion with negative electrostatic potential.

[0011] By using wavefunction analysis software to extract the characteristic parameters of the electrostatic potential distribution on the molecular surface of each ion based on wavefunction files, including the surface area and average value of the positive and negative electrostatic potential parts, key input parameters can be provided for the empirical formula of ion crystal density based on electrostatic potential. This supports the rapid prediction of the crystal density of potential ion combinations, improving the efficiency and guidance of the early screening.

[0012] Optionally, the preset empirical formula for ionic crystal density is: ; Where ρ is the predicted crystal density, M is the total mass of the ion pairs, Vm is the volume of isolated gas phase molecules, and A s + The surface area of ​​the positive electrostatic potential portion of the cation. For A s + The average value of the corresponding electrostatic potential, A s - It represents the surface area of ​​the negative electrostatic potential portion of the anion. For A s - The average electrostatic potential is represented by α, β, γ, and δ, which are empirical fitting parameters.

[0013] By employing a pre-defined empirical formula for ionic crystal density, which combines the total mass of ion pairs, the volume of isolated gas phase molecules, and multiple parameters such as the surface area and average value of the positive electrostatic potential of cations and the surface area and average value of the negative electrostatic potential of anions obtained from electrostatic potential analysis, a rapid and quantitative preliminary assessment of the potential crystal density formed by different combinations of cations and anions can be achieved. This allows for density-guided screening of candidate combinations before subsequent high-cost calculations, improving the resource allocation efficiency of the overall design process.

[0014] Optionally, the step of generating multiple initial crystal structure models includes: Select the target space group number and a Wyckoff position under the target space group, and determine the point group symmetry of the Wyckoff position; Determine the molecular point group symmetry of cations or anions in the target ion combination; Determine whether the molecular symmetry axis of the cation or anion in the target ion combination matches the point group symmetry of the Wyckoff position; When the molecular symmetry axis matches the point group symmetry, calculate the rotation matrix between the molecular symmetry axis and the symmetry axis at the Wyckoff position; The rotation matrix is ​​used to rotate one molecular symmetry axis of the cation or anion in the target ion combination to coincide with the symmetry axis of the Wyckoff position, and the centroid of the cation or anion in the target ion combination is placed at the Wyckoff site. Apply all the symmetry operations of the target space group to the cations or anions in the rotated target ion combination to generate an initial crystal structure model. Repeat the above steps to generate the multiple initial crystal structure models.

[0015] By employing an analytical approach that combines the target space group number with the Wyckoff position point group symmetry and the ion-molecule point group symmetry, and then, under matching conditions, aligning the molecular symmetry axis to the Wyckoff position symmetry axis based on the rotation matrix and applying all space group symmetry operations to generate an initial crystal structure model, a candidate structure compatible with the crystal space group can be systematically generated based on the principle of symmetry. This alleviates the inefficiency of completely random or enumeration searches and improves the orientation and rationality of the crystal structure generation stage.

[0016] Optionally, the steps for identifying the bonding states between ions include: Obtain the atomic coordinates and unit cell parameters in the second optimized crystal structure model; Calculate the distances between different atoms of anions and cations in the second optimized crystal structure model; The distance is compared to the sum of the van der Waals radii of the corresponding atoms; When the distance is less than the product of the sum of the van der Waals radii and a preset coefficient, it is determined that there is abnormal bonding.

[0017] By obtaining the atomic coordinates and unit cell parameters of the second optimized crystal structure model and calculating the atomic distance between cations and anions, then comparing them with the sum of the van der Waals radii of the corresponding atoms, and determining whether there is abnormal bonding based on a preset coefficient threshold, the chemical rationality of the optimized structure can be automatically detected. This helps to identify candidate structures that may become unstable due to unexpected strong interactions, thereby improving the chemical stability of the final screening results.

[0018] Optionally, the steps for detecting protonation include: Obtain the distance between hydrogen atoms and acceptor atoms in the second optimized crystal structure model; When the distance is less than a preset bond length threshold, protonation is determined to have occurred.

[0019] By obtaining the distance between hydrogen atoms and acceptor atoms in the second optimized crystal structure model and comparing it with a preset bond length threshold to determine whether protonation has occurred, the potential proton transfer or unexpected chemical changes in the optimized structure can be automatically identified. This helps to eliminate candidate designs that are distorted due to changes in ionic form and improves the structural chemical integrity of the final performance evaluation stage.

[0020] As can be seen from the above technical solutions, this application provides a high-throughput design method for crystal structures based on molecular or ionic configurations. The method includes: acquiring candidate ionic structure data; performing quantum chemical calculations on the ions in the candidate ionic structure data to obtain the geometric structure and wave function of the ions; performing electrostatic potential analysis on the ions according to the wave function to obtain the molecular surface electrostatic potential distribution parameters of the ions; calculating the predicted crystal density of ionic crystals formed by different combinations of cations and anions based on the molecular surface electrostatic potential distribution parameters and a preset empirical formula for ionic crystal density; selecting target ionic combinations with predicted crystal densities greater than the preset density threshold from the different combinations of cations and anions according to a preset density threshold; and matching the molecular symmetry axes of the cations and anions in the target ionic combinations with the symmetry axes of the selected Wyckoff positions in the target space group based on space group theory and Wyckoff position symmetry analysis, and applying the full symmetry of the target space group. Symmetric operations are performed to generate multiple initial crystal structure models. A machine learning potential model is used to pre-optimize these initial models, resulting in a first optimized crystal structure model. A first structure evaluation is performed on the first optimized crystal structure model to identify interionic bonding states and obtain first crystal energy and density data. Based on the results of the first structure evaluation, a first candidate crystal structure model that meets a first preset energy and density index is selected from the first optimized crystal structure models. A first-principles calculation method is then used to refine the first candidate crystal structure model, resulting in a second optimized crystal structure model. A second structure evaluation is performed on the second optimized crystal structure model to identify interionic bonding states, detect protonation phenomena, and obtain crystal energy and density data. Based on the results of the second structure evaluation, a crystal structure that meets preset energy and density performance indicators is selected from the second optimized crystal structure model, thus solving the problem of difficulty in efficiently screening and generating high-density metastable crystal structures. Attached Figure Description

[0021] To more clearly illustrate the technical solution of this application, the drawings used in the embodiments will be briefly introduced below. Obviously, for those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0022] Figure 1 A flowchart illustrating the high-throughput crystal structure design method based on molecular or ionic configuration provided in this application embodiment; Figure 2 A schematic diagram illustrating the generation of an initial crystal structure model in the high-throughput crystal structure design method based on molecular or ionic configuration provided in the embodiments of this application; Figure 3This is a schematic diagram illustrating the structure optimization and structure evaluation in the high-throughput crystal structure design method based on molecular or ionic configuration provided in the embodiments of this application. Detailed Implementation

[0023] The embodiments will now be described in detail, examples of which are illustrated in the accompanying drawings. When the following description relates to the drawings, unless otherwise indicated, the same numbers in different drawings represent the same or similar elements. The embodiments described below do not represent all embodiments consistent with this application. They are merely examples of systems and methods consistent with some aspects of this application.

[0024] To address the challenge of efficiently screening and generating high-density metastable crystal structures, see [reference needed]. Figure 1 This application provides a high-throughput design method for crystal structures based on molecular or ionic configurations, the method comprising: S100: Obtain candidate ion structure data, which includes the structural information of cations and anions.

[0025] The structural information of cations encompasses key parameters such as their chemical composition, atomic bonding, spatial configuration, and charge number. These parameters directly influence their interaction mode and binding energy with anions. Similarly, the structural information of anions includes their chemical composition, atomic arrangement, stereoconfiguration, and charge properties, such as the presence of high-energy groups like nitro and azide groups, and their spatial distribution. Systematic collection and processing of this structural information provides comprehensive and accurate foundational data for subsequent ion matching assessments and crystal structure predictions, ensuring accurate understanding of ion characteristics during high-throughput screening.

[0026] S200: Performs quantum chemical calculations on ions in candidate ion structure data to obtain the ion's geometry and wavefunction.

[0027] It should be understood that the "geometric structure" here specifically refers to the stable three-dimensional spatial arrangement of ions after optimization, including the precise coordinates of each atom, bond lengths, bond angles, and dihedral angles, etc. These parameters collectively determine the three-dimensional configuration of the ion; while the "wave function" is a mathematical function describing the state of electron motion in the ion. Through the wave function, key electronic structure information such as electron density distribution, molecular orbital energy and occupancy, and charge layout analysis can be further derived. This information is crucial for a deeper understanding of the nature of interactions between ions. Using the Gaussian 16 software package at a specific theoretical level (DFT / B3LYP / 6-31G)... The calculations under these conditions not only ensure that the obtained geometric structure is a stable point on the potential energy surface (confirmed by frequency calculation that there is no imaginary frequency), but also obtain a high-precision wave function, providing a solid quantum chemical data foundation for subsequent ion matching, crystal structure prediction and related property evaluation.

[0028] In some embodiments, the step of performing quantum chemical calculations on ions in candidate ion structure data includes: Within the framework of density functional theory, hybrid functionals and double zeta-polarized basis sets are used to perform geometric optimization calculations and frequency calculations for all ions, and the corresponding wavefunction files are output.

[0029] Among them, the "density functional theory framework" is a widely used theoretical method in current quantum chemical calculations. It uses electron density as the basic variable to describe the ground state properties of the system and approximates the Schrödinger equation by constructing exchange-correlation functionals. It can achieve a good balance between computational accuracy and computational cost, and is particularly suitable for studying the electronic structure of medium-sized and larger molecular or ionic systems. "Hybrid functionals" specifically refer to functional forms that mix a certain proportion of Hartree-Fock exchange energy and density functional theory exchange energy in the exchange energy functional. For example, in this embodiment, the B3LYP functional can be selected. This functional combines the Becke three-parameter exchange functional and the Lee-Yang-Parr correlation functional, and shows good universality and accuracy in handling molecular geometry optimization, vibrational frequencies, and energy calculations. A "double zeta valence shell polarization basis set" means that each valence shell atomic orbital in the basis functions is represented by a linear combination of two different basis functions (i.e., double zeta functions), and the basis set includes polarization functions describing the polarization effect of atomic orbitals (such as d-orbital polarization functions or p-orbital polarization functions, depending on the type of atom), for example, 6-31G. The basis set, where "6-31G" indicates that inner-shell electrons are described by a contracted Gaussian function, and valence electrons are described by a contracted function composed of a linear combination of three Gaussian functions and a diffuse function composed of a Gaussian function. "This indicates the addition of d-orbital polarization functions to heavy atoms (such as elements in the third period and above). This basis set can more accurately describe the distribution of valence electrons and the bonding between atoms. "Geometric structure optimization calculation" refers to iteratively adjusting the spatial coordinates of each atom in the ion to minimize the total energy of the system, i.e., finding the minimum point on the potential energy surface. The optimization process typically uses algorithms such as gradient descent to continuously correct the atomic coordinates until the root mean square value of the force is less than a set threshold. "Frequency calculation" involves performing simple harmonic vibration frequency analysis on the optimized structure after geometric structure optimization, calculating the second derivative matrix of the system (He... The vibrational frequencies are obtained by using the Gaussian matrix. If all calculated vibrational frequencies are positive (without imaginary frequencies), it can be confirmed that the structure is a stable point on the potential energy surface (i.e., a stable ground state configuration), rather than a transition state or an unstable saddle point. "Outputting the corresponding wavefunction file" means that after the calculation is completed, the software will generate a file containing the system's wavefunction information (such as a .fchk or .wfn format file in Gaussian software). These files store key data such as molecular orbital coefficients, occupied numbers, electron density, and layout analysis, serving as a direct data source for subsequent studies such as ion matching, interaction energy calculation, and crystal structure prediction.

[0030] S300: Perform electrostatic potential analysis on ions based on the wave function to obtain the electrostatic potential distribution parameters of the ion's molecular surface.

[0031] In some embodiments, electrostatic potential analysis includes: calculating and extracting characteristic parameters of the molecular surface electrostatic potential distribution of each ion based on a wavefunction file using wavefunction analysis software. The characteristic parameters include the surface area and average value of the positive electrostatic potential portion of the ion, and the surface area and average value of the negative electrostatic potential portion of the ion.

[0032] Specifically, the "molecular surface" here refers to the surface through which electron density isosurfaces (usually set to 0.001 electrons / Bohr) are measured. 3The van der Waals surface of a molecule, as defined by [reference needed], directly reflects the electrostatic interaction regions exhibited by the molecule. The "surface area with positive electrostatic potential" refers to the total area occupied by regions with electrostatic potential values ​​greater than zero on the aforementioned molecular surface. Its "average value" is the result of arithmetically averaging the electrostatic potential values ​​of these positive electrostatic potential regions. This parameter combination can quantitatively characterize the overall ability and intensity distribution of cations to provide positive charge interaction sites. Similarly, the "surface area with negative electrostatic potential" is the total area of ​​regions with electrostatic potential values ​​less than zero on the molecular surface. Its "average value" is the arithmetic mean of the electrostatic potential values ​​of these negative electrostatic potential regions. This set of parameters can be used to measure the comprehensive ability and intensity level of anions to provide negative charge interaction sites. By simultaneously extracting these two types of characteristic parameters for both cations and anions, the matching degree of electrostatic attraction between cations and anions can be systematically compared. For example, when the ratio of the positive electrostatic potential surface area of ​​cations to the negative electrostatic potential surface area of ​​anions is appropriate, and the absolute values ​​of their respective average values ​​are close, it is easier to form a stable crystal structure with strong electrostatic interaction and tight packing.

[0033] S400: Based on the molecular surface electrostatic potential distribution parameters and the preset empirical formula for ionic crystal density, calculate the predicted crystal density of ionic crystals formed by different combinations of cations and anions.

[0034] In some embodiments, the preset empirical formula for ionic crystal density is: ; Where ρ is the predicted crystal density; M is the total mass of the ion pairs; Vm is the volume of isolated gas phase molecules; A s + The surface area of ​​the positive electrostatic potential portion of the cation; For A s + The average value corresponding to the electrostatic potential; A s - It is the surface area of ​​the negative electrostatic potential portion of the anion; For A s - The average value of the corresponding electrostatic potential; α, β, γ, and δ are all empirical fitting parameters.

[0035] S500: Based on a preset density threshold, select target ion combinations from different combinations of cations and anions whose predicted crystal density is greater than the density threshold.

[0036] "Based on a preset density threshold" specifically refers to the fact that before performing rapid prediction of crystal packing density, researchers pre-set a specific crystal density standard based on the application requirements and performance expectations of the target energetic material, such as "greater than 1.80 g / cm³" mentioned in the example in the text. 3This threshold is not set arbitrarily, but rather takes into account the balance of energy output, safety, stability, and other aspects of energetic materials. It is usually determined based on the density range of existing high-performance energetic materials and the performance targets for new materials.

[0037] "Screening from different cation and anion combinations" refers to the various possible cation-anion pairings with potential development value obtained in the early stages of the method of this invention through quantum chemical calculations and electrostatic potential analysis. The number of these combinations can be extremely large, and directly performing subsequent complex crystal structure predictions and first-principles optimizations on all combinations would consume enormous computational resources and time. Therefore, screening by a preset density threshold involves selecting from these numerous combinations those whose crystal density, as predicted by empirical formulas, exceeds that threshold; these are the "target ion combinations." The core purpose of this step is to efficiently narrow the research scope, ensuring that subsequent, more accurate but costly computational resources are concentrated on the most promising candidate systems.

[0038] S600: Based on space group theory and Wyckoff position symmetry analysis, the molecular symmetry axes of cations and anions in the target ion combination are matched with the symmetry axes of selected Wyckoff positions in the target space group, and all symmetry operations of the target space group are applied to generate multiple initial crystal structure models.

[0039] It should be understood that, firstly, the target space group is usually based on a prediction of the type of symmetry that the target crystal may possess, or by selecting a space group commonly found in a specific type of compound as a candidate. After selecting the space group, suitable Wyckoff positions need to be chosen from it. Wyckoff positions are equivalent point systems with specific symmetries in the space group. Each Wyckoff position has its own definite symmetry characteristics, including the symmetry operations that restrict the atoms or molecules at that position.

[0040] Next, the molecular symmetry axes are matched with the Wyckoff position symmetry axes. Here, "molecular symmetry axes" refer to the symmetry axes possessed by the cations and anions in the target ion combination as independent molecules, such as the order (e.g., 2nd, 3rd, etc.) and orientation of the symmetry axis. "Wyckoff position symmetry axes," on the other hand, refer to the symmetry axis properties exhibited by the Wyckoff position under the symmetry operations of the target space group. The matching process ensures that the molecular symmetry axes of the cations and anions are compatible with the symmetry axes of the Wyckoff positions they will occupy in terms of type and orientation, allowing the molecules to "naturally" embed themselves into the lattice environment defined by the Wyckoff position without conflicting with the overall symmetry of the space group.

[0041] Then, after symmetry axis matching is completed and the cations and anions in the target ion combination are initially placed at the selected Wyckoff positions, all symmetry operations of the target space group are applied. These symmetry operations include translation, rotation, reflection, inversion, and combinations thereof. By applying these symmetry operations to the initially placed cations and anions, the positions of all atoms or molecules in the entire lattice are uniquely determined, thereby systematically generating multiple initial crystal structure models that satisfy the symmetry requirements of the target space group. These initial models form the basis for subsequent structure optimization and performance calculations.

[0042] In some embodiments, see Figure 2 The steps for generating multiple initial crystal structure models include: S610: Select the target space group number and a Wyckoff position under the target space group, and determine the point group symmetry of the Wyckoff position.

[0043] It should be understood that the Wyckoff position is an important concept in crystallography describing the equivalent position of an atom or ion in a crystal structure. It is defined by the symmetry operations of the space group. Atoms or ions belonging to the same Wyckoff position can be interconverted under crystal symmetry transformations, thus possessing identical chemical and geometric environments. Each Wyckoff position has a specific sign, multiplicity (i.e., the number of times the equivalent position recurs in a unit cell), and point group symmetry (the set of symmetry operations possessed by that position in its local environment). The specific target space group number provided above is used to select a Wyckoff position within the target space group and analyze its point group symmetry to determine the point group symmetry of the Wyckoff position.

[0044] S620: Determine the molecular point group symmetry of cations or anions in the target ion combination.

[0045] Specifically, by analyzing the molecular point group symmetry of cations or anions in the target ion combination and combining it with the point group symmetry of the Wyckoff position, all symmetry axes of the molecule are determined.

[0046] S630: Determine whether the molecular symmetry axis of the cation or anion in the target ion combination matches the point group symmetry of the Wyckoff position.

[0047] It should be understood that "matching" specifically means that the set of symmetry axes of a molecule must be a subset of the set of symmetry axes accommodated by the point group symmetry at Wyckoff positions. In other words, all symmetry axes of a molecule must be compatible with the point group symmetry at Wyckoff positions; that is, the point group symmetry at Wyckoff positions must completely encompass the symmetry axis requirements of the molecule. If a molecule has a symmetry axis that cannot be accommodated by the point group symmetry at Wyckoff positions, then it is considered a mismatch.

[0048] S640: When the molecular symmetry axis matches the point group symmetry, calculate the rotation matrix between the molecular symmetry axis and the symmetry axis at the Wyckoff position.

[0049] S650: Use a rotation matrix to rotate one molecular symmetry axis of the cation or anion in the target ion combination to coincide with the symmetry axis of the Wyckoff position, and place the centroid of the cation or anion in the target ion combination at the Wyckoff site.

[0050] S660: Apply all symmetry operations of the target space group to the cations or anions in the rotated target ion combination to generate an initial crystal structure model.

[0051] It should be understood that the target space group is... Figure 2 In the context of crystal space groups, all symmetry operations of the target space group are not merely applied once to a single cation or anion molecule in a rotational combination placed at a specific Wyckoff site. Instead, every symmetry operation within the space group (e.g., translation, rotation, reflection, inversion, and their combinations) is applied comprehensively to the initially placed molecule. Theoretically, through this series of operations, multiple molecules associated with the initial molecule are "replicated" throughout the entire periodic lattice of the crystal. These molecules generated through symmetry operations, together with the initially placed molecule, constitute the preliminary arrangement of the cation or anion in the crystal structure, forming the so-called "initial crystal structure model." This model initially demonstrates the periodic arrangement of the ion in the crystal under the selected space group and Wyckoff position conditions.

[0052] S670: Repeat the above steps to generate multiple initial crystal structure models.

[0053] It should be understood that "repeating the above steps" specifically includes performing the entire process at different Wyckoff positions within the same target space group, as well as selecting different target space groups and choosing different Wyckoff positions within each selected space group. In this way, different combinations of space groups and Wyckoff positions are systematically traversed, thereby generating multiple initial crystal structure models with different symmetries and molecular arrangements.

[0054] S710: Use a machine learning potential model to perform pre-optimization on multiple initial crystal structure models to obtain the first optimized crystal structure model.

[0055] S720: Perform a first structure evaluation on the first optimized crystal structure model to identify the bonding states between ions and obtain the first crystal energy and density data.

[0056] S730: Based on the results of the first structure evaluation, select the first candidate crystal structure model that meets the first preset energy and density index from the first optimized crystal structure model, and perform fine optimization on the first candidate crystal structure model using the first principle calculation method to obtain the second optimized crystal structure model.

[0057] Specifically, the "first optimized crystal structure model" refers to the crystal structure with relatively low energy and a reasonable structure obtained after pre-optimization using the aforementioned machine learning potential model. In the fine optimization stage, the VASP first-principles calculation method is employed, selecting exchange-correlation functionals suitable for energetic material systems (such as the PBE-D3 functional considering van der Waals interactions) and high-precision basis sets (such as PAW pseudopotentials with sufficiently high cutoff energies) to perform more rigorous and precise geometric optimization on these pre-optimized structures. This process not only further relaxes atomic positions but may also optimize lattice parameters according to research needs to fully consider the fine interactions between atoms (including covalent bonds, ionic bonds, hydrogen bonds, and van der Waals forces), ultimately obtaining the "second optimized crystal structure model" with the lowest energy and the most stable structure. This step is a precise correction of the machine learning potential pre-optimization results, ensuring that subsequent property calculations based on this crystal structure model (such as detonation performance and sensitivity parameters) have high reliability and prediction accuracy, providing a solid structural foundation for high-throughput screening and design of energetic materials. Specifically, the first preset energy and density indices are as follows: in the pre-optimized first optimized crystal structure model, the first preset energy index is set to have a formation energy lower than -0.5 eV / atom, and the lattice energy is reduced by at least 10% compared to the initial model; the first preset density index is a crystal density of not less than 1.85 g / cm³. 3This is because a lower formation energy indicates higher thermodynamic stability of the crystal structure; a significant decrease in lattice energy means more complete structural relaxation and a more equilibrium atomic arrangement; and 1.85 g / cm³ 3 The density threshold, based on the initial density screening and combined with the accuracy expectation calculated using first-principles calculations, further enhances the energy performance potential of candidate structures, laying the foundation for accurate evaluation of key performance parameters such as detonation parameters. This combination of indicators effectively screens out first-line candidate crystal structure models that possess both good stability and high density.

[0058] S800: Performs a second structure evaluation on the second optimized crystal structure model to identify interionic bonding states, detect protonation, and obtain crystal energy and density data. See details below. Figure 3 .

[0059] In some embodiments, the step of identifying the bonding state between ions includes: Obtain the atomic coordinates and unit cell parameters in the second optimized crystal structure model.

[0060] It should be understood that the "second optimized crystal structure model" specifically refers to the crystal structure model corresponding to the CONTCAR file output after structural optimization by VASP software. The core of this step is ensuring that the acquired data is an accurate result optimized through first-principles calculations. Its atomic coordinates are fractional or Cartesian coordinates within a periodic unit cell, while the unit cell parameters include lattice constants (a, b, c) and lattice angles (α, β, γ). These parameters collectively define the periodic repeating units of the crystal and are the basis for subsequent calculations of interatomic distances. The system directly reads this key information from the CONTCAR file through automated scripts, avoiding errors that may arise from manual input and ensuring data accuracy and the efficiency of the acquisition process.

[0061] Calculate the distances between different atoms of anions and cations in the second optimized crystal structure model.

[0062] It should be understood that "calculating the distances between different atoms of cations and anions in the second optimized crystal structure model" is not simply calculating the distances between all atomic pairs. Instead, it first requires identifying the types of atoms and determining their ionic properties within the crystal structure. The system distinguishes cations and anions based on the element type and its common valence states in the compound. For example, in typical oxides, metals are usually cations, and oxygen is an anion. Subsequently, for the identified cations and anions, the system iterates through all possible atomic combinations, calculating the straight-line distance between their centers of mass. This calculation needs to consider the periodicity of the unit cell; that is, when an atom is located at the unit cell boundary, it is "mirrored" into adjacent unit cells using periodic boundary conditions to accurately calculate the true distances between nearest-neighbor atoms, avoiding distance calculation errors caused by unit cell truncation.

[0063] The distance is compared to the sum of the van der Waals radii of the corresponding atoms.

[0064] It should be understood that the "van der Waals radius of the corresponding atom" refers to the van der Waals radius value of the specific atom species constituting the cation-anion pair. These radius values ​​are usually derived from authoritative chemical databases, such as the standard van der Waals radius data provided in the CRC Handbook of Chemistry and Physics. For ions, although their actual radii may vary due to the chemical bonding environment, the sum of the van der Waals radii serves as an important reference threshold for initially determining whether strong interactions exist between atoms beyond conventional ionic bonds. For example, the sum of the van der Waals radii for oxygen and hydrogen atoms is approximately 1.40 Å + 1.20 Å = 2.60 Å. If the calculated distance of an OH atom pair is much smaller than this value, it may indicate an abnormal bonding situation.

[0065] When the distance is less than the product of the sum of van der Waals radii and a preset coefficient, an abnormal bond is determined to exist.

[0066] It should be understood that the "preset coefficient" is an empirical parameter determined after extensive testing and verification, typically ranging from 0.85 to 0.95. This coefficient is set because the sum of van der Waals radii reflects the range of weak interatomic interactions (such as van der Waals forces). When the interatomic distance is significantly smaller than this range, it implies the possible existence of stronger interactions, such as covalent bonds. By introducing the preset coefficient, the threshold for "abnormality" can be defined more precisely, avoiding misjudging normal ionic bonds or weak interactions as abnormal bonding. For example, if the preset coefficient is set to 0.9, the judgment threshold for the aforementioned OH atom pair would be 2.60 Å × 0.9 = 2.34 Å. When the calculated OH distance is less than 2.34 Å, the system determines that abnormal bonding exists between the cation and anion, and this candidate structure will be flagged and excluded. This preset coefficient can be fine-tuned according to the characteristics of specific compound systems to adapt to the bonding characteristics of different types of materials, ensuring both accuracy and flexibility in judgment.

[0067] In some embodiments, the step of detecting protonation includes: Obtain the distance between hydrogen atoms and acceptor atoms in the second optimized crystal structure model.

[0068] When the distance is less than the preset bond length threshold, protonation is determined to have occurred.

[0069] It should be understood that the "second optimized crystal structure model" specifically refers to the crystal structure model presented in the CONTCAR file after VASP first-principles optimization. This model contains the precise coordinate information of all atoms in the system. In ionic crystal systems, the "acceptor atom" usually refers to the atom capable of accepting protons (H atoms). +The bond length threshold is an anion or neutral molecule, such as oxygen (as in oxides, hydroxides, and oxyanions) or nitrogen, which have lone pairs of electrons and can form stable chemical bonds with protons. The "preset bond length threshold" is not a fixed value, but rather a comprehensive setting based on the typical bond length range of covalent bonds formed by hydrogen atoms with different acceptor atoms, combined with the criteria for judging bonding interactions in computational chemistry. For example, for the OH bond, its typical covalent bond length is usually around 0.95 Å to 1.05 Å. Therefore, the system sets the OH bond length threshold at a value slightly above the upper limit of this range (e.g., 1.2 Å). When the distance between a hydrogen atom and an oxygen atom in the second optimized crystal structure model is detected to be less than this threshold, it is determined that the hydrogen atom has undergone protonation to the oxygen atom. Similarly, for the NH bond, its bond length threshold is adjusted accordingly based on the typical bond length range of the NH covalent bond (e.g., 1.00 Å to 1.10 Å). The purpose of setting this threshold is to accurately distinguish between real protonated covalent bonds and hydrogen bonds or van der Waals interactions that may exist in the crystal, ensuring the accuracy of protonation phenomenon judgment, and thus effectively identifying invalid structures that do not conform to the original ion stoichiometry or stoichiometry due to proton transfer.

[0070] S900: Based on the results of the second structure evaluation, select crystal structures that meet the preset energy and density performance indicators from the second optimized crystal structure model.

[0071] It should be understood that the "results of the second structure assessment" specifically encompass two key aspects: First, the aforementioned chemical rationality verification results, which accurately determine whether there are invalid structures in the second optimized crystal structure model due to proton transfer or other reasons that do not conform to the original ion ratio or stoichiometry, thus ensuring the chemical stability and rationality of the structure, by setting bond length thresholds (such as OH bonds, NH bonds, etc.). Second, the results of the simultaneous comprehensive performance assessment of energy and density. The system automatically extracts the final energy and crystal density data after VASP optimization and constructs a comprehensive scoring system based on preset energy threshold ranges (e.g., requiring the crystal formation energy to be lower than a certain value to ensure its thermodynamic stability or energy performance advantage) and density threshold ranges (e.g., requiring the crystal density to be higher than a certain value to meet specific material application requirements). Therefore, "screening crystal structures that meet the preset energy and density performance indicators" means that in the second optimized crystal structure model, only candidate structures that not only pass all chemical stability checks (i.e., the chemical rationality portion of the second structure evaluation results is qualified), but whose final energy value after VASP optimization falls within the preset energy performance indicator range and whose crystal density value also meets the preset density performance indicator requirements (i.e., the combined energy and density performance portion of the second structure evaluation results is also qualified) will be ultimately selected by the system and output as successful design schemes. This process ensures that the final crystal structure possesses both chemical rationality and stability, as well as energy and density characteristics that meet the design goals.

[0072] As can be seen from the above technical solutions, the embodiments of this application provide a high-throughput design method for crystal structures based on molecular or ionic configurations. The method includes: acquiring candidate ionic structure data; performing quantum chemical calculations on the ions in the candidate ionic structure data to obtain the geometric structure and wave function of the ions; performing electrostatic potential analysis on the ions according to the wave function to obtain the molecular surface electrostatic potential distribution parameters of the ions; calculating the predicted crystal density of ionic crystals formed by different combinations of cations and anions based on the molecular surface electrostatic potential distribution parameters and a preset empirical formula for ionic crystal density; selecting target ionic combinations with predicted crystal densities greater than the density threshold from different combinations of cations and anions according to a preset density threshold; and matching the molecular symmetry axes of the cations and anions in the target ionic combinations with the symmetry axes of the selected Wyckoff positions in the target space group based on space group theory and Wyckoff position symmetry analysis, and applying all the symmetry axes of the target space group. The process involves generating multiple initial crystal structure models; pre-optimizing these models using a machine learning potential model to obtain a first optimized crystal structure model; performing a first structure evaluation on the first optimized crystal structure model to identify interionic bonding states and obtain first crystal energy and density data; based on the results of the first structure evaluation, selecting first candidate crystal structure models that meet first preset energy and density indices from the first optimized crystal structure models, and performing fine optimization on the first candidate crystal structure models using first-principles calculations to obtain a second optimized crystal structure model; performing a second structure evaluation on the second optimized crystal structure model to identify interionic bonding states, detect protonation phenomena, and obtain crystal energy and density data; and based on the results of the second structure evaluation, selecting crystal structures that meet preset energy and density performance indices from the second optimized crystal structure model to address the problem of efficiently screening and generating high-density metastable crystal structures.

Claims

1. A high-throughput design method for crystal structures based on molecular or ionic configurations, characterized in that, The method includes: Acquire candidate ion structure data, which includes structural information of cations and anions; Quantum chemical calculations are performed on the ions in the candidate ion structure data to obtain the ion geometry and wavefunction; Based on the wave function, an electrostatic potential analysis is performed on the ion to obtain the molecular surface electrostatic potential distribution parameters of the ion; Based on the molecular surface electrostatic potential distribution parameters and the preset empirical formula for ionic crystal density, the predicted crystal density of ionic crystals formed by different combinations of cations and anions is calculated. Based on a preset density threshold, target ion combinations with a predicted crystal density greater than the density threshold are selected from the different combinations of cations and anions. Based on space group theory and Wyckoff position symmetry analysis, the molecular symmetry axes of the cations and anions in the target ion combination are matched with the symmetry axes of the selected Wyckoff positions in the target space group, and all symmetry operations of the target space group are applied to generate multiple initial crystal structure models. The multiple initial crystal structure models were pre-optimized using a machine learning potential model to obtain the first optimized crystal structure model. A first structure evaluation is performed on the first optimized crystal structure model to identify the bonding states between ions and obtain the energy and density data of the first crystal. Based on the results of the first structure evaluation, a first candidate crystal structure model that meets the first preset energy and density index is selected from the first optimized crystal structure model, and the first candidate crystal structure model is finely optimized using the first principle calculation method to obtain the second optimized crystal structure model. A second structure evaluation is performed on the second optimized crystal structure model to identify interionic bonding states, detect protonation, and obtain crystal energy and density data. Based on the results of the second structure evaluation, crystal structures that meet the preset energy and density performance indicators are selected from the second optimized crystal structure model.

2. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The steps for performing quantum chemical calculations on the ions in the candidate ion structure data include: Within the framework of density functional theory, hybrid functionals and double zeta-polarized basis sets are used to perform geometric optimization calculations and frequency calculations for all ions, and the corresponding wavefunction files are output.

3. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The electrostatic potential analysis includes: Based on the wavefunction file, the characteristic parameters of the molecular surface electrostatic potential distribution of each ion are calculated and extracted using wavefunction analysis software. The characteristic parameters include the surface area and average value of the part of the ion with positive electrostatic potential, and the surface area and average value of the part of the ion with negative electrostatic potential.

4. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The preset empirical formula for the density of ionic crystals is: ； Where ρ is the predicted crystal density, M is the total mass of the ion pairs, Vm is the volume of isolated gas phase molecules, and A s + The surface area of ​​the positive electrostatic potential portion of the cation. For A s + The average value of the corresponding electrostatic potential, A s - It represents the surface area of ​​the negative electrostatic potential portion of the anion. For A s - The average electrostatic potential is represented by α, β, γ, and δ, which are empirical fitting parameters.

5. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The steps for generating multiple initial crystal structure models include: Select the target space group number and a Wyckoff position under the target space group, and determine the point group symmetry of the Wyckoff position; Determine the molecular point group symmetry of cations or anions in the target ion combination; Determine whether the molecular symmetry axis of the cation or anion in the target ion combination matches the point group symmetry of the Wyckoff position; When the molecular symmetry axis matches the point group symmetry, calculate the rotation matrix between the molecular symmetry axis and the symmetry axis at the Wyckoff position; The rotation matrix is ​​used to rotate one molecular symmetry axis of the cation or anion in the target ion combination to coincide with the symmetry axis of the Wyckoff position, and the centroid of the cation or anion in the target ion combination is placed at the Wyckoff site. Apply all the symmetry operations of the target space group to the cations or anions in the rotated target ion combination to generate an initial crystal structure model. Repeat the above steps to generate the multiple initial crystal structure models.

6. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The steps for identifying the bonding state between ions include: Obtain the atomic coordinates and unit cell parameters in the second optimized crystal structure model; Calculate the distances between different atoms of anions and cations in the second optimized crystal structure model; The distance is compared to the sum of the van der Waals radii of the corresponding atoms; When the distance is less than the product of the sum of the van der Waals radii and a preset coefficient, it is determined that there is abnormal bonding.

7. The high-throughput design method for crystal structures based on molecular or ionic configurations according to claim 1, characterized in that, The steps for detecting protonation include: Obtain the distance between hydrogen atoms and acceptor atoms in the second optimized crystal structure model; When the distance is less than a preset bond length threshold, protonation is determined to have occurred.

Images