Method for constructing reactive force field based on multi-pole moment model

Abstract

The invention belongs to the field of molecular simulation, and particularly relates to a method for constructing a reactive force field based on a multi-pole moment model. Comprising the following contents: a charge penetration effect, a polarization effect, a charge transfer effect, Van der Waals interaction and a three-body potential energy form suitable for a chemical reaction; a dynamic electric multi-pole moment method suitable for chemical reaction; the invention relates to a method for constructing a potential energy surface data set with six degrees of freedom. The invention provides a simple and effective force field model for researching chemical reaction and a parameter optimization method for users. According to the method, a potential energy surface data set composed of a large number of conformations with six degrees of freedom can be generated, and a data basis is provided for model development and parameter optimization. The method is mainly used for researching the chemical reaction process and the mechanism thereof, and has a wide application prospect in the research fields of chemistry, material science and the like.

Description

TECHNICAL FIELD

[0001] The present application belongs to the field of molecular simulation, and particularly relates to a method for constructing a reactivity force field based on a multipole moment model. BACKGROUND

[0002] Molecular simulation is a method for researching chemical, material and biological systems by using computers and numerical calculation methods, has atomic spatial resolution and femtosecond time resolution, can obtain the dynamic behavior of the system in a near real environment, intuitively shows the micro mechanism, explains the experimental phenomena, and predicts the theoretical results, thereby providing strong support for disease diagnosis and treatment, innovative drug development, reaction mechanism research and chemical catalyst development.

[0003] Common molecular simulation methods include quantum chemistry method, molecular dynamics simulation and Monte Carlo simulation. The quantum chemistry method is based on quantum mechanics principle, obtains the electronic structure and energy of the system by solving the Schrodinger equation, has high calculation accuracy, can be used for researching chemical reactions, but has large calculation amount, poor scalability, and is difficult to apply to large systems (<1000 atoms). The molecular dynamics simulation is based on Newtonian mechanics, has simple form, can be applied to large systems, but the simulation process keeps the chemical bond unchanged, and cannot be used for researching chemical reactions. The Monte Carlo simulation can describe the macroscopic process (state-state process) of chemical reactions through dynamic constants, but is based on random process, cannot obtain the "real" physical process, and must obtain all the dynamic constants of reactions in advance.

[0004] The reactivity force field does not preset fixed chemical bond connection relationship, but dynamically determines the bonding relationship based on chemical environment, can describe the generation and rupture of chemical bonds during chemical reactions. Meanwhile, the reactivity force field is based on Newtonian mechanics, usually adopts the all-atom model, has small calculation amount, good scalability, and can realize the dynamic simulation of large-scale chemical reaction systems (<100,000 atoms), thereby making up for the deficiencies of the quantum chemistry method, the molecular dynamics simulation and the Monte Carlo simulation. However, the current various reactivity force field methods have deficiencies, one of which is that complex corrections based on chemical prior knowledge are adopted and there are a large number of parameters without physical meaning, and the other of which is that the point charge model lacks effective description of the complex electrostatic characteristics of molecules and the complex electrostatic interaction between molecules, which greatly limits the further development. SUMMARY

[0005] The application aims to develop a method for constructing a reactive force field based on a multipole moment model. The model uses an electric multipole moment atomic model, considers the effects of atomic charge penetration, polarization and charge transfer, and can more accurately describe the electrostatic properties of molecules and complex intermolecular electrostatic interactions, and is closer to the real chemical reaction process. At the same time, the parameters can be trained based on a large number of conformational potential energy surface data sets to obtain more general parameters. The application provides a simple and effective simulation method for users to study chemical reactions.

[0006] The technical scheme adopted by the application to achieve the above-mentioned purpose is:

[0007] A method for constructing a reactive force field based on a multipole moment model, comprising the following steps:

[0008] 1) generating a series of conformations;

[0009] 2) performing quantum chemistry calculation on the generated conformations to obtain quantum chemistry data;

[0010] 3) optimizing atomic electric multipole moment parameters based on conformation coordinates and quantum chemistry data;

[0011] 4) optimizing energy parameters based on the optimized atomic electric multipole moment parameters;

[0012] 5) performing chemical reaction dynamics simulation based on the optimized parameter values.

[0013] The step 1) is specifically:

[0014] A large number of conformational potential energy surface data sets with six degrees of freedom are generated, wherein the six degrees of freedom are six degrees of freedom of relative motion of two bodies in kinematics, including three translational degrees of freedom and three rotational degrees of freedom, the three translational degrees of freedom are specifically represented by spherical coordinates , and the three rotational degrees of freedom are specifically represented by rotation angles (α, β, γ) around the coordinate axes;

[0015] The generation of the conformation includes two ways, wherein:

[0016] Way one: a compound is split into two fragments by removing a specified chemical bond, a large number of conformations are generated by changing the relative positions of the two fragments, and the relative position changes use all six degrees of freedom

[0017] Way two: a compound is split into three fragments by removing two specified chemical bonds, a large number of conformations are generated by changing the relative positions of the three fragments, and the relative position changes use three translational degrees of freedom for each chemical bond The subscripts 1 and 2 of the variables represent whether the variable corresponds to the first chemical bond or the second chemical bond.

[0018] The step 3) is specifically:

[0019] The atom regenerates the atomic electric multipole moment according to the atomic environment in which the atom is located, wherein the electric monopole moment q i of the atom i is:

[0020]

[0021] Wherein, represents the charge flow, and represents that there is a corresponding charge amount flowing from the atom j to the atom i;

[0022] The charge flow is contributed by three parts, the first part is the charge flow caused by the difference in element types between atoms, that is:

[0023]

[0024] Wherein, the intermediate variable f q,i is:

[0025]

[0026] The second part is the charge flow caused by the difference in atomic environment between atoms of different element types, that is:

[0027]

[0028] The third part is the charge flow caused by the difference in atomic environment between atoms of the same element type, that is:

[0029]

[0030] Wherein, A1-A8, k1-k8, r1-r8 and B0-B3 are parameters, mp q , dp q and qp q are monopole moment, dipole moment and quadrupole moment in atomic distribution multipole moment, respectively, for distinguishing different atomic environments, r ij is the distance between atoms;

[0031] The electric dipole moment μ i of the atom i is:

[0032]

[0033] The electric quadrupole moment Q i of the atom i is:

[0034]

[0035] Wherein δ αβ is the Dirac function, δ αβ =1 when α=β; δαβ = 0, val i represents the number of valence electrons of atom i, and respectively represent the atomic distribution of generating atomic electric dipole moment and electric quadrupole moment, specifically:

[0036]

[0037] wherein A 10 , k 10 and r 10 are parameters.

[0038] The atomic distribution multipole moment of atom i is expressed as:

[0039]

[0040] wherein α, β = x, y, z represent three components of the dipole moment vector and nine components of the quadrupole tensor, represents the atomic distribution, and its expression is:

[0041]

[0042] wherein A9, k9 and r9 are parameters;

[0043] The atomic distribution dipole moment and quadrupole moment respectively adopt the modulus dp cp，i of the atomic distribution dipole moment and the Frobenius norm qp cp，i of the atomic distribution quadrupole moment as descriptors, wherein

[0044] The step 4) is specifically: calculating the expression of the reactivity force field model energy is:

[0045] E = E cp + E vdW + E pol + E ct + E tri

[0046] wherein E cp , E vdW , E pol , E ct and E tri respectively represent the charge penetration energy, the van der Waals energy, the polarization energy, the charge transfer energy and the three-body potential.

[0047] The charge penetration energy is:

[0048] E cp = ∑ i ∑ j＞i f cp，ij · Ecp，ij

[0049] where E cp，ij represents the charge penetration energy between atom i and atom j, and its expression is:

[0050]

[0051] where r ij is the interatomic distance, γ z is a coefficient, and Z is the nuclear charge number. M = [q, μ x , μ y , μ z , Q xx , Q xy , Q xz , …, Q zz ] T is the atomic electric multipole moment, and q, μ and Q represent the atomic monopole, dipole and quadrupole moments, respectively. and are the multipole interaction matrices.

[0052] f cp，ij represents the environmental factor of the charge penetration effect, and its expression is:

[0053]

[0054] where represents the atomic distribution, and its expression is:

[0055]

[0056] where the upper index C0represents the index, and C0-C5, A 11 , k 11 and r 11 are parameters, mp cp，ij , dp cp，ij and qp cp，ij are the monopole, dipole and quadrupole moments in the atomic distribution multipole moment, respectively.

[0057] The atomic distribution multipole moment of atom i is expressed as:

[0058]

[0059] where α, β = x, y, z represent the three components of the dipole moment vector and the nine components of the quadrupole moment tensor, and the atomic distribution dipole moment and the atomic distribution quadrupole moment are described by the modulus dp cp，i of the atomic distribution dipole moment and the Frobenius norm qp cp，i of the atomic distribution quadrupole moment, respectively, where The combination rule is:

[0060]

[0061] where val represents the number of valence electrons.

[0062] The van der Waals energy is:

[0063]

[0064] where ε ij and σ ij represent the potential well depth and the minimum energy distance, respectively, and δ and γ are buffer constants.

[0065] The polarization energy is:

[0066]

[0067] where, represents the multipole-multipole interaction matrix of the polarization effect, represents the induced dipole moment of the polarization effect, which is expressed as:

[0068]

[0069] where represents the polarizability.

[0070] The charge transfer energy is:

[0071]

[0072] where, represents the multipole-multipole interaction matrix of the charge transfer effect, represents the induced dipole moment of the charge transfer effect, which is expressed as:

[0073]

[0074] where, represents the charge transfer rate.

[0075] The three-body potential is:

[0076] E tri =∑ i ∑ j＞i ∑ k＞j E ijk

[0077] where E ijk represents the interaction between atom i, atom j, and atom k, and its specific form is:

[0078] E ijk =f3(r ij ,rik ,θ jik )+f3(r ij ,r jk ,θ ijk )+f3(r ik ,r jk ,θ jki )

[0079]

[0080] where A θ and θ0 represent the force constant and equilibrium angle of angular vibration, k b and r0 represent the force constant and equilibrium distance of bond dissociation.

[0081] The present application has the following advantages and benefits:

[0082] 1. The reactive force field model has no any artificial correction based on chemical prior knowledge, and has a physical image more consistent with the real chemical reaction process. The total energy includes the charge penetration effect, charge transfer effect and polarization effect based on the atomic model of electric multipole moment, as well as the van der Waals interaction and three-body potential. The charge penetration effect represents the phenomenon of electron cloud overlap between atoms, which provides the bonding energy between atoms; the van der Waals interaction provides the repulsive energy between atoms at close range; the charge transfer effect describes the charge transfer phenomenon between atoms; and the polarization effect describes the electron cloud deformation under the external complex electrostatic environment. The present application makes corrections to the above energy forms based on the chemical environment, so that it can be applied to the chemical reaction process.

[0083] 2. The reactive force field model generates atomic electric multipole moments based on the atomic environment before simulating the energy and force of each step, wherein the atomic point charge is obtained by the charge flow based on the atomic distribution multipole moment, the atomic dipole moment is obtained by the vector sum of the atomic distribution, and the atomic quadrupole moment is obtained by the tensor sum of the atomic distribution. The dynamic electric multipole moment method developed by the present application makes the electric multipole moment atomic model applicable to the chemical reaction process.

[0084] 3. The present application can generate a potential energy surface dataset composed of a large number of conformations with six degrees of freedom, providing a data basis for model development and parameter optimization. BRIEF DESCRIPTION OF DRAWINGS

[0085] Figure 1 is a specific flowchart for implementing the present application;

[0086] Figure 2 is a schematic diagram of six degrees of freedom in the case of removing one chemical bond;

[0087] Figure 3 is a schematic diagram of six degrees of freedom in the case of removing two chemical bonds. DETAILED DESCRIPTION

[0088] The application will be described in further detail below with reference to the drawings and embodiments.

[0089] As shown in the drawings, a method for constructing a reactive force field based on a multipole moment model includes the following steps: Figure 1

[0090] (1) generating a series of conformations by a conformation generation program;

[0091] (2) performing quantum chemical calculation on the conformations obtained in step (1) and extracting quantum chemical calculation results by a data processing program;

[0092] (3) inputting the conformation coordinates of step (1) and the quantum chemical data of step (2) into a parameter optimization program to optimize atomic electric multipole moment parameters;

[0093] (4) fixing the parameters obtained in step (3) and optimizing energy parameters by a parameter optimization program;

[0094] (5) realizing chemical reaction kinetics simulation.

[0095] In step (1), the conformation generation program can generate a potential energy surface dataset composed of a large number of conformations with six degrees of freedom; the six degrees of freedom adopt six degrees of freedom of relative motion of two bodies in kinematics, including three translational degrees of freedom and three rotational degrees of freedom, the three translational degrees of freedom are specifically represented by spherical coordinates , and the three rotational degrees of freedom are specifically represented by rotation angles (α, β, γ) around the coordinate axes; the conformation generation program supports two conformation generation modes: the first mode splits a compound into two fragments by removing a specified chemical bond, and then generates a large number of conformations by changing the relative positions of the two fragments, wherein the relative position changes use all six degrees of freedom The second mode splits a compound into three fragments by removing two specified chemical bonds, and then generates a large number of conformations by changing the relative positions of the three fragments, wherein the relative position changes use three translational degrees of freedom for each chemical bond The six degrees of freedom of the above two modes can be coupled, so that for the first mode, a total of conformations can be generated; for the second mode, a total of conformations can be generated, wherein the six factors respectively represent the number of points on the respective degrees of freedom.

[0096] In step (4), the energy expression of the reactive force field model is:

[0097] E = E cp + E vdW + E pol + E ct + E​tri (1)

[0098] where E cp , E vdW , E pol , E ct and E tri represent the charge penetration energy, the van der Waals energy, the polarization energy, the charge transfer energy and the three-body potential, respectively;

[0099] The specific expression of the charge penetration energy is:

[0100] E cp =∑ i ∑ j＞i f cp，ij ·E cp，ij (2)

[0101] where E cp，ij represents the charge penetration energy between atom i and atom j, and its specific expression is:

[0102]

[0103] where r ij is the interatomic distance; γ z is a coefficient to enhance the charge penetration energy to be applicable to the chemical reaction system; Z is the nuclear charge number; M = [q, μ x , μ y , μ z , Q xx , Q xy , Q xz , …, Q zz ] T is the atomic electric multipole moment, and q, μ and Q represent the atomic monopole moment, the dipole moment and the quadrupole moment, respectively; and are the multipole moment interaction matrices for the charge penetration effect between the atomic nucleus and the electron, and between the electrons, respectively;

[0104] f cp，ij represents the environmental factor of the charge penetration effect, and is used to describe different chemical environments, and its specific expression is:

[0105]

[0106] where represents the atomic distribution, and its specific expression is:

[0107]

[0108] In the formula, C0-C5, A 11 , k 11 and r 11 are parameters; mpcp，ij , dp cp，ij and qp cp，ij are the monopole, dipole and quadrupole moments in the atom- distributed multipole expansion, respectively, to distinguish different chemical environments;

[0109] The atom-distributed multipole moments of atom i are expressed as:

[0110]

[0111] where α, β = x, y, z represent the three components of the dipole vector and the nine components of the quadrupole tensor; the atom-distributed dipole and quadrupole moments are respectively represented by the modulus of the atom-distributed dipole and the Frobenius norm of the atom-distributed quadrupole as descriptors to maintain the consistency of the atomic environment under different coordinate systems; the combination rule is:

[0112]

[0113] where val represents the number of valence electrons;

[0114] The specific expression of the van der Waals energy is:

[0115]

[0116] where ε ij and σ ij represent the potential well depth and the minimum energy distance, respectively, and δ and γ are buffer constants;

[0117] The specific expression of the polarization energy is:

[0118]

[0119] where represents the multipole interaction matrix of the polarization effect, represents the induced dipole moment of the polarization effect, and its specific expression is:

[0120]

[0121] where represents the polarizability;

[0122] The specific expression of the charge transfer energy is:

[0123]

[0124] where represents the multipole interaction matrix of the charge transfer effect, The induced dipole moment representing the charge transfer effect, which is specifically expressed as:

[0125]

[0126] wherein represents the charge transfer rate;

[0127] The three-body potential is specifically expressed as:

[0128] E tri =∑ i ∑ j＞i ∑ k＞j E ijk (13)

[0129] wherein

[0130] E ijk = f3(r ij ,r ik ,θ jik ) + f3(r ij ,r jk ,θ ijk ) + f3(r ik ,r jk ,θ jki ) (14) and

[0131]

[0132] wherein A θ and θ0 respectively represent the force constant and the equilibrium angle of the angular vibration, k b and R0 represent the force constant and the equilibrium distance of the bond dissociation;

[0133] In step (3), the method for generating atomic electric multipole moments of the reactive force field model is a dynamic atomic environment-based electric multipole moment method, specifically, before simulating the energy and force of each step, the atomic electric multipole moment is regenerated according to the atomic environment of the atom;

[0134] The electric monopole moment q i of atom i is specifically expressed as:

[0135]

[0136] wherein represents the charge flow, representing the corresponding charge amount flowing from atom j to atom i; the charge flow is contributed from three parts, the first part of the contribution is the charge flow caused by the difference in element types between atoms:

[0137]

[0138] wherein

[0139]

[0140] The second contribution is due to the charge flow between atoms of different element types caused by the difference in atomic environment:

[0141]

[0142] The third contribution is due to the charge flow between atoms of the same element type caused by the difference in atomic environment:

[0143]

[0144] where A1-A8, k1-k8, r1-r8 and B0-B3 are parameters, mp q , dp q and qp q are the monopole, dipole and quadrupole moments of the atom distribution multipole, respectively, to distinguish different atomic environments;

[0145] The atom distribution multipole of atom i is expressed as:

[0146]

[0147] where α,β=x,y,z represent the three components of the dipole vector and the nine components of the quadrupole tensor, represents the atom distribution, and its specific expression is:

[0148]

[0149] where A9, k9 and r9 are parameters; similarly, the atom distribution dipole and quadrupole use the modulus and the Frobenius norm of the atom distribution quadrupole as descriptors to maintain the consistency of the atomic environment under different coordinate systems;

[0150] The electric dipole moment μ i of atom i is specifically expressed as:

[0151]

[0152] The electric quadrupole moment Q i of atom i is specifically expressed as:

[0153]

[0154] where val i represents the number of valence electrons of atom i, and The atomic distribution of generating atomic dipole moment and atomic quadrupole moment, both of which adopt the same form in this model, are expressed as:

[0155]

[0156] where A 10 , k 10 and r 10 are parameters.

[0157] Example 1

[0158] This embodiment takes C2H6 molecule as an example to illustrate the construction of potential energy surface dataset composed of a large number of different C2H6 conformations, and then optimize the parameters of the reactive force field based on the dataset, and finally realize the reaction dynamics simulation. The specific steps are as follows:

[0159] (1) Generating conformations

[0160] First, draw the C2H6 molecule using chemical graphics software (such as GaussView software, etc.), and save its structure information as an xyz format file.

[0161] Then, set the configuration file of the conformation generation program, the main contents of which include the path of the xyz file, the method of quantum chemical calculation basis set and the conformation generation method. Among them, the conformation generation method adopts a self-defined format, and for the case of removing a chemical bond, the specific setting method is as follows:

[0162] eulermodenameeul__1bondatom1 atom2

[0163]

[0164] Where euler and eul_1bond are fixed keywords, modename can be freely set to distinguish different dissociation channels, and atom1 and atom2 are the two atoms of the broken bond. The next 18 parameters correspond to 6 degrees of freedom in groups of 3 (such as Figure 2 ). Take the translational degree of freedom R as an example, its corresponding R0R d R N , where R0 represents the initial value, R d represents the interval of values, and R N represents the number of points, then the value of the degree of freedom R of the i-th conformation is R i =R0+R d *i, i=0→R N In addition, the file name of the conformation coordinate file will contain its serial number in six degrees of freedom, and the actual value of each degree of freedom of the current conformation can be derived from the file name, which is convenient for screening and analysis. The above six degrees of freedom are coupled, so a total of one chemical bond, it is set as:

[0165] spheric modename sph__2bond

[0166]

[0167]

[0168] where spheric and sph_2bond are fixed keywords, atom1 and atom2 are two atoms of the first chemical bond to be broken, and the following 9 numbers correspond to 3 degrees of freedom respectively, atom3 and atom4 are two atoms of the second chemical bond to be broken, and the following 9 numbers correspond to 3 degrees of freedom respectively (e.g. Figure 3 ), the values of which are consistent with the removal of one chemical bond. The above six degrees of freedom are coupled, and a total of conformations can be generated.

[0169] Finally, the above xyz coordinate file (reference structure) and configuration file are input files, and a series of conformations are generated by the conformation generation program. The coordinates of the generated conformations are stored in an xyz format coordinate file and a com format coordinate file, one conformation per file, where the former is the input file format for the subsequent parameter optimization program, and the latter is the input file format for the subsequent quantum chemical calculation. The program also generates Gaussian software running commands with com files as input, all commands are written into one file, which is convenient for running all quantum chemical calculations with simple command lines. In addition, the program supports excluding conformations with minimum atomic distance less than a certain value to save computing resources, and these conformations will usually lead to errors or non-convergence in quantum calculation. The Gaussian software is a commonly used calculation software in the field of quantum chemistry research.

[0170] (2) Generate data

[0171] First, run the command file obtained in step (1), and calculate the quantum chemical data of all conformations by the Gaussian program in order.

[0172] Then, set the configuration file of the data processing program.

[0173] Finally, the above configuration file is input file, and the Gaussian program output results are arranged by the data processing program.

[0174] The program will extract energy information from the Gaussian output file and store the energy information and the corresponding coordinate file name in the energy file mol_names_all_energy (mol represents the molecule name, which can be set in the configuration file) as an index file for reading coordinates in subsequent parameter optimization. The energy is the relative energy relative to the reference structure, and the reference structure file name is identified by reference. At the same time, the program will extract other data from the Gaussian output file and store these data and the corresponding coordinate file name in the data file mol_data_all (mol represents the molecule name, consistent with the energy file). The current program supports extracting atomic Mulliken charge, molecular multipole moment and atomic force. In addition, the program will also recombine a series of xyz format coordinate files generated by the conformation generation program into one file and store it in the coordinate file mol_txyz_all (mol represents the molecule name, consistent with the energy file), to avoid the extremely time-consuming reading of a large number of small files by the subsequent parameter optimization program. The files generated by the conformation generation program and the data processing program need to be written in a certain format and should be avoided to be manually modified.

[0175] (3) Optimize parameters

[0176] First, set the configuration file of the parameter optimization program, which mainly includes the paths of the energy file, coordinate file, data file and parameter file, as well as the energy interval, optimization step number, optimization method and optimization content. The energy interval setting program only reads the conformations whose energy relative to the reference structure is within the energy interval. The optimization method currently supported by the program is the simulated annealing algorithm. The optimization content currently supported by the program is the optimization of atomic electric multipole moment and the optimization of energy. The parameter file can freely select the parameters to be optimized and limit the parameter range.

[0177] Then, optimize the atomic electric multipole moment parameters with Mulliken charge as the target.

[0178] Finally, fix the atomic electric multipole moment parameters and optimize the energy parameters with energy as the target. The parameter optimization program uses GPU acceleration to realize parameter optimization based on a large number of conformations.

[0179] (4) Dynamics simulation

[0180] Based on the parameters in step (3), the reaction dynamics simulation is realized.

Claims

1. A method for constructing a reactive force field based on a multipole moment model, characterized in that, Includes the following steps: 1) Generate a series of conformations; 2) Perform quantum chemical calculations on the generated conformation to obtain quantum chemical data; 3) Optimize atomic electric multipole moment parameters based on conformational coordinates and quantum chemical data; 4) Based on the optimized atomic electric multipole moment parameters, the energy parameters are optimized; 5) Perform chemical reaction kinetics simulations based on the optimized parameter values.

2. The method for constructing a reactive force field based on a multipole moment model according to claim 1, characterized in that, Step 1) specifically refers to: A potential energy surface dataset consisting of numerous conformations with six degrees of freedom is generated. These six degrees of freedom are derived from the kinematics of relative motion between two bodies, including three translational degrees of freedom and three rotational degrees of freedom. Specifically, the three translational degrees of freedom are expressed in spherical coordinates. The three rotational degrees of freedom are specifically represented by rotation angles (α, β, γ) around the coordinate axes; The generation of the conformation includes two methods, wherein: Method 1: By removing a specified chemical bond, the compound is broken down into two fragments. Changing the relative positions of the two fragments generates a large number of conformations, where the relative position changes utilize all six degrees of freedom. Method 2: By removing two specified chemical bonds, the compound is broken down into three fragments. Changing the relative positions of the three fragments generates a large number of conformations, where the relative position changes utilize three translational degrees of freedom for each chemical bond. The subscripts 1 and 2 indicate whether the variable corresponds to the first or the second chemical bond.

3. The method for constructing a reactive force field based on a multipole moment model according to claim 1, characterized in that, Step 3) specifically refers to: Atoms regenerate their atomic electric multipole moments based on their atomic environment, where the electric singlet moment q of atom i is... i for: in, This indicates the flow of charge, representing the flow of a corresponding amount of charge from atom j to atom i. The charge flow originates from three contributions. The first contribution is due to the charge flow caused by the difference in elemental types between atoms, namely: Among them, the intermediate variable f q,i for: The second contribution is due to the charge flow caused by the differences in atomic environments between atoms of different element types, namely: The third contribution is due to the charge flow caused by differences in atomic environments between atoms of the same element type, namely: Where A1-A8, k1-k8, r1-r8, and B0-B3 are all parameters, mp q dp q and qp q These are the singlet moment, dipole moment, and quadrupole moment in the atomic distribution multipole moment, used to distinguish different atomic environments. r ij The distance between atoms; The electric dipole moment μ of atom i i for: The electric quadrupole moment Q of atom i i for: Where δ αβ Let δ be the Dirac function, and when α = β, δ αβ =1; when α≠β, δ αβ =0,val i This represents the number of valence electrons in atom i. and The atomic distributions representing the generating atomic electric dipole moments and electric quadrupole moments are as follows: Among them, A 10 k 10 and r 10 For parameters.

4. The method for constructing a reactive force field based on a multipole moment model according to claim 3, characterized in that, The atomic distribution multipole moment of atom i is expressed as: Where α,β=x,y,z represent the three components of the dipole moment vector and the nine components of the quadrupole moment tensor. The atomic distribution is expressed as follows: Where A9, k9, and r9 are parameters; The atomic distributed dipole moment and the quadrupole moment are respectively expressed using the modulus dp of the atomic distributed dipole moment. cp,i The Frobenius norm qp of the atomic distribution quadrupole moment cp,i As a descriptor, where, 5. The method for constructing a reactive force field based on a multipole moment model according to claim 1, characterized in that, Step 4) specifically involves calculating the energy of the reactive force field model using the following expression: E=E cp +E vdW +E pol +E ct +E tri Among them, E cp E vdW E pol E ct and E tri These represent charge penetration energy, van der Waals energy, polarization energy, charge transfer energy, and three-body potential, respectively.

6. The method for constructing a reactive force field based on a multipole moment model according to claim 5, characterized in that, The charge penetration energy is: AND cp =∑ i ∑ j＞i f cp,ij ·AND cp,ij Among them, E cp,ij The charge penetration energy between atoms i and j is expressed as follows: Where, r ij γ is the interatomic distance. z Where Z is the nuclear charge number, and M = [q, μ] is the coefficient. x ,μ y ,μ z Q xx Q xy Q xz ,…,Q zz ] T The atomic electric multipole moment is represented by q, μ, and Q, which represent the atomic monopole moment, dipole moment, and tetrapole moment, respectively. and The multipole moment interaction matrix; f cp,ij The environmental factor representing the charge penetration effect is expressed as follows: in, The atomic distribution is expressed as follows: Wherein, the superscript C0 represents the exponent, C0-C5, A 11 k 11 and r 11 All are parameters, mp cp,ij dp cp,ij and qp cp,ij These are the singlet moment, dipole moment, and quadrupole moment in the multipole moments of atomic distribution; The atomic distribution multipole moment of atom i is expressed as: Where α,β=x,y,z represent the three components of the dipole moment vector and the nine components of the quadrupole moment tensor, respectively. The atomic distributed dipole moments and quadrupole moments are represented by the modulus dp of the atomic distributed dipole moments. cp,i The Frobenius norm qp of the atomic distribution quadrupole moment cp,i As a descriptor, where, The combination rules are as follows: Where val represents the number of valence electrons.

7. The method for constructing a reactive force field based on a multipole moment model according to claim 5, characterized in that, The van der Waals energy is: Where, ε ij and σ ij δ and γ represent the potential well depth and minimum energy distance, respectively, and are buffer constants.

8. A method for constructing a reactive force field based on a multipole moment model according to claim 5, characterized in that, The polarization energy is: in, The multipole moment interaction matrix representing the polarization effect. The induced dipole moment representing the polarization effect is expressed as: in This represents the polarizability.

9. A method for constructing a reactive force field based on a multipole moment model according to claim 5, characterized in that, The charge transfer energy is: in, The multipole moment interaction matrix represents the charge transfer effect. The induced dipole moment representing the charge transfer effect is expressed as follows: in, This represents the charge transfer rate.

10. A method for constructing a reactive force field based on a multipole moment model according to claim 5, characterized in that, The three-body potential is: AND tri =∑ i ∑ j＞i ∑ k＞j AND ijk Among them, E ijk The interaction between atoms i, j, and k is represented in the following form: E ijk =f3(r ij ,r ik ,θ jik )+f3(r ij ,r jk ,θ ijk )+f3(r ik ,r jk ,θ jki ) Among them, A θ θ0 and k represent the force constant and equilibrium angle of angular vibration, respectively. b r0 represents the force constant and equilibrium distance for bond dissociation.

Images