JP2021032617A - Crystal structure calculation method, crystal structure calculation program and crystal structure calculation device - Google Patents

Abstract

To provide a crystal structure calculation method capable of precisely calculating the crystal structure of a target substance, a crystal structure calculation program, and a crystal structure calculation device.SOLUTION: The crystal structure calculation method includes the steps of: sequentially changing an initial crystal structure of a target substance, which includes a lattice vector variable representing a reciprocal lattice vector of a unit lattice and an atomic position variable representing a position vector of each unit cell in the unit lattice of each atom, within a range of Cartesian product set of lattice vector variable and atomic position variable of target substance; calculating a free energy of a pair of modified lattice vector variable and atomic position variable; and in the continuous function of free energy whose domain is a continuous subset satisfying formula (1) among the Cartesian products, calculating the pair of lattice vector variable and the atomic position variable of the target substance based on whether or not the calculated free energy is the minimum.SELECTED DRAWING: Figure 1

Description

The present invention relates to a crystal structure calculation method, a crystal structure calculation program, and a crystal structure calculation device, and in particular, a crystal structure calculation method, a crystal structure calculation program, and a crystal that calculate the crystal structure of a target substance based on lattice energy and free energy. Regarding the structure calculation device.

The crystal structure of a substance is one of the main factors that determine the properties of the substance. For example, the mobility of organic semiconductors and the solubility of pharmaceuticals largely depend on the crystal structure of the substance. The crystal structure of this substance is generally determined by single crystal X-ray diffraction measurement, and it has been difficult to predict its crystal structure before synthesizing a new substance. If it becomes possible to predict the crystal structure on a computer before synthesizing a new substance, the development of the new substance can be greatly streamlined.

Therefore, as a technique for predicting a crystal structure on a computer, for example, Non-Patent Document 1 proposes a lattice energy minimization (LEM) method for searching for a crystal structure that minimizes lattice energy. .. This LEM method is to search for a crystal structure in which _{the total potential energy of atoms (lattice energy Ust} ) is the minimum or minimum, after ignoring the kinetic energy of the atoms by assuming that the atoms are stationary. is there.

Report on the sixth blend test of organic crystal structure prediction methods, Anthony M. et al. Really et al. , Acta Cryst. (2016). B72, P.I. 439-P. 459

However, since the LEM method ignores the thermal vibration and zero-point vibration of atoms, it is not possible to reproduce thermal phenomena such as thermal expansion and structural transition due to temperature, and there is a difference from the actual crystal structure. May occur. In addition, when a substance that can obtain only a few crystal structures (polymorphs) in an actual experiment is simulated by the LEM method, hundreds or more crystal structures (polymorphs) may be obtained. There was also a problem that a crystal structure far apart was calculated.

The present invention has been made to solve such a conventional problem, and provides a crystal structure calculation method, a crystal structure calculation program, and a crystal structure calculation device for calculating the crystal structure of a target substance with high accuracy. With the goal.

In the crystal structure calculation method according to the present invention, the initial crystal structure of the target substance including a lattice vector variable indicating the basic lattice vector of the unit cell and an atomic position variable indicating the position vector of each atom in the unit cell is obtained from the target substance. Sequentially change within the range of the direct product set of the lattice vector variable and the atomic position variable, calculate the free energy of the changed lattice vector variable and the atomic position variable set, and continuously satisfy the following equation (1) in the direct product set. In a continuous function of free energy whose definition area is a subset, the set of the lattice vector variable and the atomic position variable of the target substance is calculated based on whether or not the calculated free energy shows the minimum.

Here, the lattice vector variable is sequentially changed as an independent variable within the range of the direct product set, and the atomic position variable indicating the minimum lattice energy with respect to the changed lattice vector variable is calculated as the dependent variable, thereby being included in the subset. The set of lattice vector variable and atomic position variable is searched, and the free energy of the set of lattice vector variable and atomic position variable included in the subset is calculated sequentially, and the free energy shows the minimum lattice vector variable and atomic position variable. It is preferable to calculate the set of.

In addition, the lattice vector variable and the atomic position variable are sequentially changed as independent variables within the range of the Cartesian product, and the changed set of the lattice vector variable and the atomic position variable is included in the subset based on Lagrange's undetermined multiplier method. Moreover, it is possible to determine whether or not the subset has the minimum free energy.

In addition, the free energy can be calculated based on the natural frequency of phonons calculated by pseudo-harmonic oscillator approximation from quantum chemistry calculation or classical force field calculation.

In addition, an empirical formula for calculating the free energy based on the values of the lattice vector variable and the atomic position variable included in the crystal structure information obtained in the past is calculated in advance, and the free energy is calculated based on the empirical formula. You can also.

The crystal structure calculation program according to the present invention causes a computer to execute each step of the crystal structure calculation method described in any of the above.

The crystal structure calculation device according to the present invention determines the initial crystal structure of the target substance including a lattice vector variable indicating the basic lattice vector of the unit cell and an atomic position variable indicating the position vector of each atom in the unit cell. The variable change part that is sequentially changed within the range of the direct product set of the atomic position variable and the atomic position variable, and the free energy of the set of the changed lattice vector variable and the atomic position variable are calculated, and the continuous that satisfies the following equation (1) in the direct product set. Minimal condition judgment unit that calculates the set of lattice vector variable and atomic position variable of the target substance based on whether the calculated free energy shows the minimum in the continuous function of free energy whose definition area is a specific subset. It is equipped with.

According to the present invention, the reciprocal lattice vector variable and the atomic position variable are used in a continuous function of free energy whose domain is a continuous subset satisfying the following equation (1) among the direct product sets of the reciprocal lattice vector variable and the atomic position variable. Since the set of the lattice vector variable and the atomic position variable of the target substance is calculated based on whether or not the free energy calculated by sequentially changing shows the minimum, the crystal structure calculation for calculating the crystal structure of the target substance with high accuracy. It becomes possible to provide a method, a crystal structure calculation program, and a crystal structure calculation device.

It is a block diagram which shows the structure of the crystal structure calculation apparatus which concerns on Embodiment 1 of this invention. It is a flowchart which shows the crystal structure calculation method which concerns on Embodiment 1. It is a figure which shows the state of searching the set of the lattice vector variable and the atomic position variable of the target substance in Embodiment 1. It is a figure which shows the state of searching the minimum point of the lattice energy. It is a figure which shows the state of searching the set of the lattice vector variable and the atomic position variable of the target substance by the conventional LEM method. It is a figure which shows the state of calculating the crystal structure candidate. It is a figure which shows the minimum point in the continuous function of free energy. It is a block diagram which shows the structure of the crystal structure calculation apparatus which concerns on Embodiment 2. It is a flowchart which shows the crystal structure calculation method which concerns on Embodiment 2. It is a figure which shows the mode of searching the set of the lattice vector variable and the atomic position variable of the target substance in Embodiment 2. It is a block diagram which shows the main part of the structure of the crystal structure calculation apparatus which concerns on Embodiment 3. It is a block diagram which shows the main part of the structure of the crystal structure calculation apparatus which concerns on Embodiment 4. FIG. It is a block diagram which shows the main part of the structure of the crystal structure calculation apparatus which concerns on Embodiment 5. It is a graph which shows the result of having calculated the change of the lattice constant with respect to temperature by using the crystal structure calculation method of Embodiment 1. It is a graph which shows the result of having calculated the change of the free energy with respect to the lattice constant by using the crystal structure calculation method of Embodiment 1.

Hereinafter, embodiments of the present invention will be described with reference to the accompanying drawings.
Embodiment 1
FIG. 1 shows the configuration of the crystal structure calculation device 1 according to the first embodiment of the present invention. The crystal structure calculation device 1 has an initial structure generation unit 2, and the initial structure generation unit 2 includes a lattice energy calculation unit 3, a minimum condition determination unit 4, a free energy calculation unit 5, a minimum condition determination unit 6, and an output unit 7. Are connected in sequence. Further, the minimum condition determination unit 4 is connected to the atomic position variable change unit 8, the minimum condition determination unit 6 is connected to the lattice vector variable change unit 9, and the atomic position variable change unit 8 and the lattice vector variable change unit 9 are respectively. It is connected to the lattice energy calculation unit 3. Then, the variable changing unit 10 is formed by the atomic position variable changing unit 8 and the lattice vector variable changing unit 9.
Further, the initial structure generation unit 2, the lattice energy calculation unit 3, the minimum condition determination unit 4, the free energy calculation unit 5, the minimum condition determination unit 6, the atomic position variable change unit 8, and the lattice vector variable change unit 9 have a control unit 11. It is connected, and the operation unit 12 and the storage unit 13 are connected to the control unit 11. Further, the display unit 14 is connected to the output unit 7 of the crystal structure calculation device 1.

The initial structure generation unit 2 includes a lattice vector variable indicating the reciprocal lattice vector of the unit cell and an atomic position variable indicating the position vector of each atom in the unit cell based on the molecular structure information of the target substance input from the user. Generates the initial crystal structure of the target substance containing.
Here, the molecular structure information is information on the molecular structure of the target substance, and examples thereof include the type of atoms contained in the molecule of the target substance, the number of atoms, and the type of chemical bond between atoms.

In the initial crystal structure of the target substance generated by the initial structure generation unit 2, the atomic position variable changing unit 8 fixes the lattice vector variable within the range of the reciprocal lattice set of the reciprocal lattice vector variable and the atomic position variable of the target substance, and atomizes the atom. Change the position variable sequentially. Here, the Cartesian product shows all combinations of the lattice vector variable and the atomic position variable of the target substance, and can be defined based on the molecular structure information.

The lattice energy calculation unit 3 is an atomic position variable set by the initial structure generation unit 2, the atomic position variable change unit 8 or the lattice vector variable change unit 9 based on the lattice energy calculation information of the target substance input from the user. And the lattice vector variable set is used to calculate the lattice energy.
Here, the lattice energy calculation information is information for calculating the lattice energy, and includes, for example, molecular structure information and force field parameters, and may further include environmental pressure (pressure of the environment in which crystals are formed). preferable. As the force field parameter, for example, when the covalent bond energy in the lattice energy is expressed by f (r) = A (r−r ₀ ) ² (A is the spring constant, r is the distance between atoms, and r ₀ is. _{A and r 0} of the distance between atoms in the equilibrium state) can be mentioned.

The minimum condition determination unit 4 determines whether or not the lattice energy calculated by the lattice energy calculation unit 3 indicates the minimum in the continuous function of the lattice energy calculated by sequentially changing the atomic position variables while fixing the lattice vector variables. To judge. As a result, among the direct product sets of the lattice vector variable and the atomic position variable of the target material, the atomic position variable change unit 8 is added to the subset in which the set of the lattice vector variable and the atomic position variable satisfying the following equation (1) is continuously continued. It is determined whether or not the set of the lattice vector variable and the atomic position variable set in is included.

Here, _Ust is the lattice energy of the crystal, the vectors a, b and c are the lattice vector variables, the vector r _i (i = 1, 2, ..., N) is the atomic position variable, and N is the atom in the unit cell. Indicates the number of.

The lattice vector variable changing unit 9 sequentially changes the value of the lattice vector variable within the range of the direct product set for the set of the lattice vector variable and the atomic position variable determined by the minimum condition determination unit 4 to indicate the minimum lattice energy.
In this way, the lattice vector variable changing unit 9 sequentially changes the lattice vector variable as an independent variable within the range of the direct product set of the lattice vector variable and the atomic position variable of the target substance, and the minimum condition determination unit 4 is changed. By calculating the atomic position variable whose lattice energy is the minimum with respect to the lattice vector variable as the dependent variable, the pair of the lattice vector variable and the atomic position variable included in the subset is searched.

Based on the free energy calculation information of the target substance input from the user, the free energy calculation unit 5 is determined by the minimum condition determination unit 4 that the lattice energy indicates the minimum, that is, it is determined that the lattice energy is included in the subset. The free energy is calculated from the set of the lattice vector variable and the atomic position variable.
Here, the free energy calculation information is information for calculating the free energy, and includes, for example, molecular structure information, force field parameters, environmental temperature (temperature of the environment in which crystals are formed), and environmental pressure (environmental pressure). It is preferable to include the pressure of the environment in which crystals are formed).

In the minimum condition determination unit 6, whether or not the free energy sequentially calculated by the free energy calculation unit 5 indicates the minimum in the continuous function of the free energy calculated from the set of the lattice vector variable and the atomic position variable included in the subset. Based on the above, the set of the lattice vector variable and the atomic position variable indicating the crystal structure of the target substance is calculated.

The output unit 7 outputs the crystal structure of the target substance determined by the minimum condition determination unit 6 to the display unit 14.
The control unit 11 controls each unit of the crystal structure calculation device 1 based on an instruction input from the operation unit 12 by the user.

The operation unit 12 is for inputting a command from the user, and can be formed of a button, a touch panel, a keyboard, a mouse, a trackball, or the like.
The storage unit 13 stores an operation program or the like, and stores a hard disk, a flexible disk, MO, MT, RAM, CD-ROM, DVD-ROM, SD card, CF card, a recording medium such as a USB memory, a server, or the like. Can be used.

The initial structure generation unit 2, the lattice energy calculation unit 3, the minimum condition determination unit 4, the free energy calculation unit 5, the minimum condition determination unit 6, the atomic position variable change unit 8, the lattice vector variable change unit 9, and the variable change unit 10. Consists of a CPU and an operation program for causing the CPU to perform various processes, but these may be configured by a digital circuit.
The display unit 14 includes a display device and displays the crystal structure of the target substance output from the output unit 7.

Next, the crystal structure calculation method of the first embodiment will be described with reference to the flowchart shown in FIG. The crystal structure calculation program of the present invention is a program that causes the computer of the crystal structure calculation device 1 to execute the following crystal structure calculation method.

First, when the molecular structure information of the target substance is input to the initial structure generation unit 2, the initial structure generation unit 2 sets a lattice vector variable and an atomic position variable based on the molecular structure information in step S1. To generate the initial crystal structure of the target substance. In general, the crystal structure consists of a lattice vector variable (with 6 degrees of freedom) indicating the x, y, and z components of the three reciprocal lattice vectors of the unit cell, and x, y, of the position vector of each atom in the unit cell. It can be represented by an atomic position variable (degree of freedom is 3N) indicating the z component, and the initial structure generation unit 2 generates the initial crystal structure of the target substance by selecting a set of the lattice vector variable and the atomic position variable. be able to.
Although the lattice vector variable excludes the degree of freedom of rotation, it may include the degree of freedom of rotation. Further, in the atomic position variable, N indicates the number of atoms in the unit cell. Moreover, the molecular structure information may be saved in advance.

Here, the initial structure generation unit 2 selects a set of the lattice vector variable and the atomic position variable within the range of the Cartesian product set of the lattice vector variable and the atomic vector variable of the target substance calculated based on the molecular structure information.
For example, as shown in FIG. 3, the initial structure generation unit 2 may randomly select a set V1 of a lattice vector variable and an atomic position variable within the range of the Cartesian product T1 of the lattice vector variable and the atomic vector variable of the target substance. it can. Here, the contour line L indicates the distribution of the lattice energy in the Cartesian product T1. The initial structure generation unit 2 outputs the set V1 of the selected lattice vector variable and the atomic position variable to the lattice energy calculation unit 3.
It is preferable that the initial structure generation unit 2 selects the set V1 of the lattice vector variable and the atomic position variable by limiting the selection range and the like based on the molecular structure information.

When the set V1 of the lattice vector variable and the atomic position variable selected by the initial structure generation unit 2 is input to the lattice energy calculation unit 3, the lattice energy calculation unit 3 performs the target substance input by the user in step S2. Based on the lattice energy calculation information of, the lattice energy is calculated using the set V1 of the lattice vector variable and the atomic position variable. The calculated lattice energy is output from the lattice energy calculation unit 3 to the minimum condition determination unit 4. Lattice energies include, for example, bond distance energy between two chemically bonded atoms, bond angle energy between three chemically bonded atoms, biplane angle energy between four chemically bonded atoms, and molecules of two non-chemically bonded atoms. Lattice energy and the like can be mentioned.
The lattice energy calculation information may be input by the user together with the molecular structure information, or may be stored in advance.

When the lattice energy calculated by the lattice energy calculation unit 3 is input to the minimum condition determination unit 4, the minimum condition determination unit 4 calculates by sequentially changing the atomic position variables while fixing the lattice vector variables in step S3. In the continuous function of the lattice energy to be performed, it is determined whether or not the lattice energy calculated by the lattice energy calculation unit 3 shows the minimum.
Here, since the lattice energy calculated from the set V1 of the lattice vector variable and the atomic position variable is the one input first, it cannot be determined whether or not to indicate the minimum, and the minimum condition determination unit 4 indicates the minimum. It is determined that there is no such determination, and the determination result is output to the atomic position variable changing unit 8.

When the determination result of the minimum condition determination unit 4 is input to the atomic position variable change unit 8, the process proceeds to step S4, and the atomic position variable change unit 8 has only atomic position variables from the set V1 of the lattice vector variable and the atomic position variable. Is changed in the range of the reciprocal lattice set T1 to obtain the set V2 of the reciprocal lattice vector variable and the atomic position variable, and the set V2 of the reciprocal lattice vector variable and the atomic position variable is output to the lattice energy calculation unit 3.

When the set V2 of the lattice vector variable and the atomic position variable generated by the atomic position variable changing unit 8 is input to the lattice energy calculation unit 3, the lattice energy calculation unit 3 returns to step S2 and returns to the lattice vector variable and the atomic position. The lattice energy is calculated using the set of variables V2.
In this way, the atomic position variable changing unit 8 sequentially changes the value of the atomic position variable while fixing the lattice vector variable, and the lattice energies of the set V2, V3 and V4 of the lattice vector variable and the atomic position variable are latticed. It is input from the energy calculation unit 3 to the minimum condition determination unit 4.

The minimum condition determination unit 4 minimizes the lattice energy calculated by the lattice energy calculation unit 3 in the continuous function of the lattice energy in step S3 based on the lattice energy so far input sequentially from the lattice energy calculation unit 3. Judge whether or not to indicate.

For example, as shown in FIG. 4, when the lattice energies U1 to U4 are sequentially input from the lattice energy calculation unit 3 so far, the minimum condition determination unit 4 sets the smallest value in the lattice energies U1 to U4. The indicated lattice energy U3 is determined as the minimum point. Here, the lattice energies U1 to U4 need not be calculated by changing the atomic position variables at regular intervals, and can be calculated by changing the atomic position variables at different intervals. For example, the lattice energy U2 can be calculated by changing the value of the atomic position variable to a large value with respect to the lattice energy U1, and the lattice energies U3 and U4 can be calculated by changing the value of the atomic position variable by a small amount.
The minimum condition determination unit 4 calculates the derivative of the lattice energy sequentially input from the lattice energy calculation unit 3, that is, the above equation (1), and determines the lattice energy U3 whose derivative indicates zero as the minimum point. You can also.

In this way, as shown in FIG. 3, the minimum condition determination unit 4 is connected to the subset T2 of the Cartesian product set T1 in which the set of the lattice vector variable and the atomic position variable satisfying the above equation (1) are continuously continued. , It is determined whether or not the set V1 to V4 of the lattice vector variable and the atomic position variable is included. The minimum condition determination unit 4 searches for the minimum point of the lattice energy, ignoring the thermal vibration and the zero-point vibration of the atom, in the same manner as the conventional lattice energy minimization (LEM) method.
As a result, when the lattice vector variable is fixed, the minimum condition determination unit 4 searches for an atomic position variable representing the center of atomic vibration, that is, a search for an atomic position variable included in a subset T2 in which free energy can be calculated. be able to.

Then, the atomic position variable is sequentially changed in the atomic position variable changing unit 8 until the minimum condition determination unit 4 determines that the lattice energy indicates the minimum. When the minimum condition determination unit 4 determines that the lattice energy of the set V3 of the lattice vector variable and the atomic position variable indicates the minimum, for example, the minimum condition determination unit 4 outputs the set V3 of the lattice vector variable and the atomic position variable to the free energy calculation unit 5.

As a result, the free energy calculation unit 5 calculates the free energy in step S5 using the set V3 of the lattice vector variable and the atomic position variable based on the free energy calculation information of the target substance input from the user. The free energy calculation information may be input by the user together with the molecular structure information, or may be stored in advance.
Here, the free energy calculation unit 5 can calculate based on Helmholtz's free energy F = U _st + U _ph −TS _ph or Gibbs free energy G = F + PV. Here, U _ph is the internal energy of phonons (atomic vibrations), T is absolute temperature, S _ph is the entropy of phonons, P is pressure, and V is volume. Based on擬調sum oscillator _{approximation, U ph,} using the natural vibration of the phonon set {[nu _i}, it is expressed as _{U ph = Σ i (hν i} / 2) coth (hν i / 2k B T) , _{S ph} is expressed as _{S ph = Σ i {-k ·} ln [2sinh (hν i / 2k B T)] + (hν i / 2T) coth (hν i / 2k B T)}.

The free energy calculation unit 5 outputs the calculated free energy to the minimum condition determination unit 6. When the free energy calculated by the free energy calculation unit 5 is input to the minimum condition determination unit 6, the minimum condition determination unit 6 is a set of the lattice vector variable and the atomic position variable included in the subset T2 in step S6. In the calculated continuous function of free energy, that is, the continuous function of free energy having the subset T2 as the domain, it is determined whether or not the free energy calculated by the free energy calculation unit 5 shows the minimum.
When the minimum condition determination unit 6 determines that the free energy calculated by the free energy calculation unit 5 does not show the minimum, the determination result is output to the lattice vector variable change unit 9.

When the determination result of the minimum condition determination unit 6 is input to the lattice vector variable change unit 9, the process proceeds to step S7, and the lattice vector variable change unit 9 changes only the lattice vector variable from the set V3 of the lattice vector variable and the atomic position variable. Is changed within the range of the reciprocal lattice set T1 to obtain a set V5 of the reciprocal lattice vector variable and the atomic position variable, and the set V5 of the reciprocal lattice vector variable and the atomic position variable is output to the lattice energy calculation unit 3.

In this way, when the set V5 of the lattice vector variable and the atomic position variable generated by the lattice vector variable changing unit 9 is input to the lattice energy calculation unit 3, the lattice energy calculation unit 3 returns to step S2 and the lattice The lattice energy is calculated using the set V5 of the vector variable and the atomic position variable. Subsequently, in the same manner as described above, the minimum condition determination unit 4 sequentially changes the atomic position variables while fixing the lattice vector variables until it is determined that the lattice energy calculated by the lattice energy calculation unit 3 indicates the minimum. ..
That is, the lattice vector variable is sequentially changed as an independent variable in the range of the Cartesian product set T1, and the atomic position variable indicating the minimum lattice energy with respect to the changed lattice vector variable is calculated as the dependent variable. As a result, the set of the lattice vector variable and the atomic position variable included in the subset T2 can be reliably calculated.

For example, when the minimum condition determination unit 4 determines that the lattice energy of the set V6 of the lattice vector variable and the atomic position variable indicates the minimum, the free energy calculation unit 5 frees the set V6 of the lattice vector variable and the atomic position variable. Calculate the energy.
The minimum condition determination unit 6 is a lattice calculated by the free energy calculation unit 5 in step S6 based on the free energies of the set V3 and V6 of the lattice vector variable and the atomic position variable sequentially input from the free energy calculation unit 5. It is determined whether or not the free energy of the set V6 of the vector variable and the atomic position variable shows the minimum in the continuous function of the free energy.

When the minimum condition determination unit 6 determines that the free energy shows the minimum in the set V6 of the lattice vector variable and the atomic position variable, for example, the set V6 of the lattice vector variable and the atomic position variable is used as the optimized crystal structure. To do.
As in the case of the minimum condition determination unit 4, the minimum condition determination unit 6 does not need to change the set of the lattice vector variable and the atomic position variable included in the subset T2 at regular intervals, and the lattice included in the subset T2. The minimum free energy can be determined by changing the pair of vector variable and atomic position variable at different intervals. Further, the minimum condition determination unit 6 can also determine the free energy in which the derivative of the free energy sequentially input from the free energy calculation unit 5 indicates zero as the minimum point.

Conventionally, in the LEM method, as shown in FIG. 5, the crystal structure of the target substance is a set of discrete lattice vector variables and atomic position variables V1', V2', and V3', which show the minimum value at the contour line L of the lattice energy. Calculate as. Here, since the LEM method considers the atom to be stationary, there is no concept of temperature, and only the crystal structure at 0K can be predicted (it is not 0K because the zero-point vibration is ignored). There is. At this time, even when the free energy is calculated for the set V1', V2'and V3'of the lattice vector variable and the atomic position variable to determine the crystal structure of the target substance, the set V1 of the lattice vector variable and the atomic position variable is determined. Since the', V2'and V3' are discrete, it is not possible to reproduce continuous structural changes such as thermal expansion, and the problem of excess polymorphism cannot be solved.

The present invention is a free energy minimization (FEM) method for searching a set V6 of a lattice vector variable and an atomic position variable showing a minimum in a continuous function of free energy. Here, according to the second law of thermodynamics, it is equivalent that the crystal is in a thermal equilibrium state under isothermal isobaric conditions and the Gibbs free energy of the crystal is minimal. Therefore, by searching for a crystal structure in which the Gibbs free energy is minimized, it is possible to calculate a crystal structure closer to reality in consideration of thermal vibration and zero-point vibration.
Specifically, the crystal structure of the target substance can be calculated with high accuracy in consideration of the magnitude of phonons (atomic vibrations) in continuous temperature changes. Further, in the prediction of the crystal structure, since the crystal structure at room temperature is important, the calculation of the crystal structure at a finite temperature instead of 0K is extremely important. In the present invention, it is possible to reproduce thermal expansion and structural transition due to temperature with high accuracy. Furthermore, in the present invention, since only a crystal structure that can stably exist at a finite temperature can be obtained, it is possible to suppress the calculation of an excessive number of crystal structures (polymorphs).

Such a method of calculating the crystal structure is performed in parallel or continuously. That is, as shown in FIG. 3, the initial structure generation unit 2 randomly selects a new set V1a of the lattice vector variable and the atomic position variable in the range of the direct product set T1, and repeats the above steps S2 to S7. , The minimum condition determination unit 6 can sequentially calculate the set V2a of the lattice vector variable and the atomic position variable indicating the minimum free energy in the subset T2 as the optimized crystal structure.
For example, as shown in FIG. 6, when 1000 initial crystal structures are randomly generated in the initial structure generation unit 2, the minimum condition determination unit 6 calculates 1000 optimized crystal structures showing the minimum free energy. Will be done.

Then, the minimum condition determination unit 6 combines 1000 optimized crystal structures having the same crystal structure into one, and calculates, for example, crystal structure candidates R1, R2 ... R9 of nine target substances. .. As a result, as shown in FIG. 7, the crystal structure candidates R1 to R9 showing the minimum points in the continuous function f of the free energy with the subset T2 as the domain can be calculated. At this time, when the number of crystal structure candidates R1 to R9 is large, the minimum condition determination unit 6 determines the calculated crystal structure candidates R1 to R9 whose free energy is equal to or less than a predetermined value, for example, the crystal structure candidate R4. , R5 and R6 are calculated as the crystal structure of the target substance.
As described above, since the minimum condition determination unit 6 excludes the crystal structure candidates R1 to R9 having a large free energy value, the crystal structures R4, R5 and R6 of the target substance can be calculated with higher accuracy.

The minimum condition determination unit 6 outputs the calculated crystal structures R4, R5 and R6 of the target substance to the display unit 14 via the output unit 7, and the crystal structures R4, R5 and R6 of the target substance are displayed on the display unit 14. To. Thereby, for example, based on the crystal structures R4, R5 and R6 of the target substance displayed on the display unit 14, the X-ray analysis data obtained in the experiment can be analyzed with high accuracy.

According to this embodiment, the free energy of the set of the lattice vector variable and the atomic position variable is calculated, and the free energy of the direct product set T1 whose domain is the continuous subset T2 satisfying the above equation (1). In the continuous function f, since it is determined whether or not the calculated free energy shows the minimum, the crystal structure of the target substance can be calculated with high accuracy.
Specifically, the lattice vector variable is sequentially changed as an independent variable in the range of the direct product set T1, and the atomic position variable indicating the minimum lattice energy with respect to the changed lattice vector variable is calculated as the dependent variable. The pair V3 and V6 of the lattice vector variable and the atomic position variable included in T2 can be surely searched. Further, in order to sequentially calculate the free energy of the set of the lattice vector variable and the atomic position variable included in the subset T2 and calculate the set V6 of the lattice vector variable and the atomic position variable whose free energy indicates the minimum, the target substance. The crystal structure of can be calculated with higher accuracy.

Embodiment 2
In the first embodiment, the minimum condition determination unit 6 sequentially calculates the free energies of the set V3 and V6 of the lattice vector variable and the atomic position variable included in the subset T2, and the lattice vector indicating that the free energy is the minimum. The set V6 of the variable and the atomic position variable was calculated, but based on whether or not the free energy of the set of the lattice vector variable and the atomic position variable shows the minimum in the continuous function f of the free energy with the subset T2 as the domain. The crystal structure of the target substance may be calculated, and the present invention is not limited to this.

For example, as shown in FIG. 8, the variable change unit 22 is arranged in place of the variable change unit 10 of the first embodiment, the minimum condition determination unit 23 is arranged in place of the minimum condition determination units 4 and 6, and the grid is further arranged. Except for the energy calculation unit 3 and the free energy calculation unit 5, the deviation direction calculation unit 21 can be newly arranged. Here, the deviation direction calculation unit 21, the variable change unit 22, and the minimum condition determination unit 23 are sequentially connected and arranged between the initial structure generation unit 2 and the output unit 7.

The shift direction calculation unit 21 calculates the shift direction for changing the set of the lattice vector variable and the atomic position variable based on the lattice energy and the free energy. Specifically, the Lagrange function L is calculated based on the Lagrange's undetermined multiplier method represented by the following equation (2) using the set of the lattice vector variable and the atomic position variable, the lattice energy calculation information, and the free energy calculation information. , The deviation direction in which the Lagrange function L becomes the smallest is calculated when the values of the lattice vector variable and the atomic position variable are slightly changed in a plurality of directions.

Here, L is the Lagrange function, G is the free energy of the target substance, the vectors a, b and c are lattice vector variables, the vector r _i (i = 1, 2, ..., N) is the atomic position variable, x _i. , Y _i , z _i are x, y, z components of the vector r _i , N is the number of atoms in the unit cell of the target substance, λ is the undetermined multiplier of Lagrange, and U _st is the lattice energy of the target substance.

The variable changing unit 22 changes the values of the lattice vector variable and the atomic position variable within the range of the direct product set T1 of the lattice vector variable and the atomic position variable based on the deviation direction calculated by the deviation direction calculation unit 21.

In the minimum condition determination unit 23, the set of the lattice vector variable and the atomic position variable changed by the variable change unit 22 is included in the subset T2 and is the minimum in the subset T2 based on the above equation (2). Determine if it has free energy. That is, the minimum condition determination unit 23 searches for a set of a lattice vector variable and an atomic position variable in which the Lagrange function L of the above equation (2) indicates a stop point (minimum point, maximum point, inflection point, etc.). Then, the minimum condition determination unit 23 calculates a set of a lattice vector variable and an atomic position variable included in the subset T2 and having the minimum free energy in the subset T2 as the crystal structure of the target substance.

Next, the crystal structure calculation method of the second embodiment will be described with reference to the flowchart shown in FIG.
First, as in the first embodiment, when the molecular structure information of the target substance is input to the initial structure generation unit 2 by the user, the initial structure generation unit 2 performs the lattice in step S1 based on the molecular structure information. Select a set of vector variables and atomic position variables to generate the initial crystal structure of the target substance. For example, as shown in FIG. 10, the initial structure generation unit 2 randomly selects a set V1 of a lattice vector variable and an atomic position variable within the range of the Cartesian product set T1 of the lattice vector variable and the atomic vector variable of the target substance. The initial structure generation unit 2 outputs the set V1 of the selected lattice vector variable and the atomic position variable to the shift direction calculation unit 21.

When the set V1 of the lattice vector variable and the atomic position variable selected by the initial structure generation unit 2 is input to the deviation direction calculation unit 21, the deviation direction calculation unit 21 is based on the above equation (2) in step S21. , Calculate the deviation direction for changing the set V1 of the lattice vector variable and the atomic position variable. Specifically, the shift direction calculation unit 21 calculates the Lagrange function L by slightly changing the value of the set V1 of the lattice vector variable and the atomic position variable in a plurality of directions, and the Lagrange function L becomes the smallest. The direction is determined to be the deviation direction N1. The deviation direction calculation unit 21 outputs the calculated deviation direction N1 to the variable change unit 22.

When the deviation direction N1 calculated by the deviation direction calculation unit 21 is input to the variable change unit 22, the variable change unit 22 sets the set V1 of the lattice vector variable and the atomic position variable in the range of the direct product set T1 in step S22. Change toward the shift direction N1. For example, the variable changing unit 22 largely changes the set V1 of the lattice vector variable and the atomic position variable to the set V2 of the lattice vector variable and the atomic position variable. The variable change unit 22 outputs the changed set V2 of the lattice vector variable and the atomic position variable to the minimum condition determination unit 23.
In this way, the variable changing unit 22 changes the set V1 of the lattice vector variable and the atomic position variable based on the deviation direction N1 calculated using the Lagrange function L, so that the lattice vector having lower lattice energy and free energy It is possible to efficiently search the set V2 of the variable and the atomic position variable.

When the set V2 of the lattice vector variable and the atomic position variable changed by the variable change unit 22 is input to the minimum condition determination unit 23, the minimum condition determination unit 23 sets the lattice vector variable and the atomic position variable in step S23. The Lagrange function L is calculated based on the above equation (2) using V2, and it is determined whether or not the Lagrange function L indicates a stop point. As a result, it is determined whether or not the set V2 of the lattice vector variable and the atomic position variable is included in the subset T2 and has the minimum free energy in the subset T2.
When the minimum condition determination unit 23 determines that the set V2 of the lattice vector variable and the atomic position variable does not indicate the stop point of the Lagrange function L, the minimum condition determination unit 23 outputs the determination result to the deviation direction calculation unit 21.

When the determination result of the minimum condition determination unit 23 is input to the deviation direction calculation unit 21, the process returns to step S21, and the deviation direction calculation unit 21 sets the lattice vector variable and the atomic position variable based on the above equation (2). The deviation direction N2 for changing V2 is calculated.
In this way, it is determined that the value of the Lagrange function L calculated by the minimum condition determination unit 23 using the set of the lattice vector variable and the atomic position variable indicates the stop point in the Lagrange function L calculated so far. Until then, the variable changing unit 22 sequentially changes the set V2 and V3 of the lattice vector variable and the atomic position variable toward the deviation directions N2 and N3. The minimum condition determination unit 23 does not need to change the set of the lattice vector variable and the atomic position variable at regular intervals, but changes the set of the lattice vector variable and the atomic position variable at different intervals to change the stop point of the Lagrange function L. It can be determined.

Then, when the value of the Lagrange function L of the set V4 of the lattice vector variable and the atomic position variable is determined to indicate a stop point in the Lagrange function L calculated so far, the minimum condition determination unit 23 determines the value. The set V4 of the lattice vector variable and the atomic position variable is used as the optimized crystal structure.
In this way, since the minimum condition determination unit 23 determines based on the Lagrange's undetermined multiplier method, a set of a lattice vector variable and an atomic position variable showing the minimum in the continuous function f of free energy with the subset T2 as the domain. V4 can be searched efficiently.

Such a method for calculating the crystal structure is performed in parallel or continuously as in the first embodiment. Then, the minimum condition determination unit 23 sequentially determines the optimized crystal structure, and the crystal structure of the target substance is extracted from the plurality of optimized crystal structures.

According to the present embodiment, the minimum condition determination unit 23 includes the set V1 to V4 of the lattice vector variable and the atomic position variable in the subset T2 and in the subset T2 based on the Lagrange's undetermined multiplier method. In order to determine whether or not it has the minimum free energy, it is possible to efficiently search the set V4 of the lattice vector variable and the atomic position variable that satisfy the condition.

Embodiment 3
In the above embodiments 1 and 2, the free energy is preferably calculated based on quantum chemical calculation or classical force field calculation.
For example, as shown in FIG. 11, the parameter storage unit 31 can be newly arranged in the first embodiment.

The parameter storage unit 31 is connected to the free energy calculation unit 5 and stores force field parameters calculated based on quantum chemical calculation or classical force field calculation. The parameter storage unit 31 can store, for example, the force field parameters of the target substance and the force field parameters of a substance similar to the target substance. Here, the force field parameters can be calculated by quantum chemistry calculation based on the molecular structure information. The force field parameters can also be calculated by classical force field calculation based on the crystal structure information obtained in the past.

As a result, when the set of the lattice vector variable and the atomic position variable included in the subset T2 is input from the minimum condition determination unit 4, the free energy calculation unit 5 forces the force field of the target substance stored in the parameter storage unit 31. Search for force field parameters of a substance similar to the parameter or target substance. Subsequently, the free energy calculation unit 5 uses the set of the force field parameter stored in the parameter storage unit 31 and the lattice vector variable and the atomic position variable input from the minimum condition determination unit 4 to approximate the phonon by the pseudo-harmonic oscillator. The natural frequency is calculated, and the free energy of the set of the lattice vector variable and the atomic position variable is calculated based on the natural frequency of the phonon.
Then, the calculated free energy is output from the free energy calculation unit 5 to the minimum condition determination unit 6, and it is determined whether or not the free energy shows the minimum in the continuous function f of the free energy with the subset T2 as the domain. Will be done.

In this way, since the parameter storage unit 31 stores the force field parameters calculated based on the quantum chemistry calculation or the classical force field calculation, the free energy calculation unit 5 sets the force field parameters stored in the parameter storage unit 31. Free energy can be calculated accurately based on this.

According to the present embodiment, since the free energy calculation unit 5 calculates the natural frequency of the phonon from the quantum chemistry calculation or the classical force field calculation based on the pseudo-harmonic oscillator approximation, it is based on the natural frequency of the phonon. The free energy can be calculated accurately.

Embodiment 4
In the third embodiment, the free energy calculation unit 5 calculates the free energy based on the natural frequency of the phonon calculated from the quantum chemistry calculation or the classical force field calculation. However, if the free energy can be calculated, Well, it's not limited to this.
For example, as shown in FIG. 12, the empirical formula storage unit 41 can be arranged in place of the parameter storage unit 31 of the third embodiment.

The empirical formula storage unit 41 stores the empirical formula for calculating the free energy based on the values of the lattice vector variable and the atomic position variable included in the crystal structure obtained in the past. Examples of the empirical formula include a formula when the value of free energy is calculated in the past from the values of the lattice vector variable and the atomic position variable for a crystal structure similar to the target substance.

As a result, when the set of the lattice vector variable and the atomic position variable included in the subset T2 is input from the minimum condition determination unit 4, the free energy calculation unit 5 is based on the empirical formula stored in the empirical formula storage unit 41. Then, the free energy of the set of the lattice vector variable and the atomic position variable is calculated.
Then, the calculated free energy is output from the free energy calculation unit 5 to the minimum condition determination unit 6, and it is determined whether or not the free energy shows the minimum in the continuous function f of the free energy with the subset T2 as the domain. Will be done.

In this way, in order for the empirical formula storage unit 41 to store the empirical formula for calculating the free energy, the free energy calculation unit 5 promptly stores the free energy based on the empirical formula stored in the empirical formula storage unit 41. Can be calculated.

According to the present embodiment, since the free energy calculation unit 5 calculates an empirical formula based on the values of the lattice vector variable and the atomic position variable included in the crystal structure information obtained in the past, it is based on the empirical formula. Free energy can be calculated quickly.

Embodiment 5
In the above-described first to fourth embodiments, when the X-ray analysis data of the target substance has already been obtained, the minimum condition determination unit refers to the X-ray analysis data and the lattice vector variable and the atom of the target substance. It is preferable to calculate the set of positional variables.
For example, as shown in FIG. 13, the X-ray data storage unit 51 can be newly arranged in the first embodiment.

The X-ray data storage unit 51 is connected to the minimum condition determination unit 6 and stores the X-ray analysis data of the target substance obtained by the experiment. Examples of the X-ray analysis data include the values of the lattice vector variable and the atomic position variable of the target substance obtained by the X-ray analysis.

As a result, when the free energy is input from the free energy calculation unit 5, the minimum condition determination unit 6 acquires the X-ray analysis data of the target substance from the X-ray data storage unit 51 and refers to the X-ray analysis data. Calculate the set of the lattice vector variable and the atomic position variable of the target substance. Specifically, the minimum condition determination unit 6 contains the values of the lattice vector variable and the atomic position variable of the free energy input from the free energy calculation unit 5 and the X-ray analysis data stored in the X-ray data storage unit 51. Calculate the degree of agreement with the values of the lattice vector variable and the atomic position variable. Subsequently, the minimum condition determination unit 6 indicates that the value obtained by adding the degree of coincidence to the free energy input from the free energy calculation unit 5 is the minimum in the continuous function f, that is, the continuous function f of the free energy to which the degree of coincidence is added. Judge whether or not. At this time, it is preferable that the minimum condition determination unit 6 adds a predetermined weight to the free energy and the degree of coincidence.

Then, the minimum condition determination unit 6 calculates the set of the free energy lattice vector variable and the atomic position variable determined to indicate the minimum as the crystal structure of the target substance. As a result, the set of the lattice vector variable and the atomic position variable of the target substance can be calculated with higher accuracy.

According to the present embodiment, the minimum condition determination unit 6 calculates the set of the lattice vector variable and the atomic position variable by referring to the X-ray analysis data as well as the free energy, so that the crystal structure of the target substance is shown. The set of the lattice vector variable and the atomic position variable can be calculated with higher accuracy. Further, the minimum condition determination unit 6 calculates the degree of coincidence between the set of the lattice vector variable and the atomic position variable changed in the range of the reciprocal lattice set T2 and the X-ray analysis data. The set of the lattice vector variable and the atomic position variable indicating the crystal structure of the material can be calculated with higher accuracy.

In the above embodiments 1 to 5, the minimum condition determination unit does not need to change the set of the lattice vector variable and the atomic position variable included in the subset T2 in order, and the set of the lattice vector variable and the atomic position variable. Can be changed in various ways to search for the minimum point of the continuous function f.

(Example 1)
The crystal structure of argon was calculated using the crystal structure calculation method of the first embodiment. As the crystal structure, the change of the lattice constant from 0K to 90K was calculated. The result is shown in FIG.
In addition, the change in free energy with respect to the lattice constant was calculated for each temperature from 0K to 90K. The result is shown in FIG.
(Comparative Example 1)
The crystal structure of argon was calculated using the LEM method. As the crystal structure, the lattice constant was calculated. The result is shown in FIG.
In addition, the change in lattice energy with respect to the lattice constant was calculated for each temperature from 0K to 90K. The result is shown in FIG.

From the results shown in FIG. 14, in Example 1, the difference in lattice constant with respect to the experimental value of thermal expansion of argon (Dobbs & Jones, 1957) was 0.7% or less, whereas in Comparative Example 1, 1.6%. It was found that a value close to the experimental value can be calculated. Further, it was found that in Comparative Example 1, the temperature axis did not exist, whereas in Example 1, the rate of thermal expansion accompanying the temperature rise was in good agreement with the experimental value. Further, in Example 1, the minimum point of free energy disappears and expands infinitely at 87 K or higher. It was found that this temperature is an intermediate value between the melting point (83.8K) and the boiling point (87.3K) of argon, and it can be reproduced that the crystal state is not stable above the boiling point. As described above, by using the FEM method, it was possible to reproduce the thermal expansion which is a thermal phenomenon, and to appropriately exclude the crystal structure which is not stable at a finite temperature.

From the results shown in FIG. 15, it was found that in Example 1, the curvature of the free energy decreased and the magnitude of the atomic vibration increased as the temperature increased. At this time, the minimum point of free energy could be calculated when the temperature was 0K to 80K, but the minimum point could not be calculated when the temperature was 90K. As a result, unstable structures at finite temperatures could be excluded, and an excessive number of crystal structures (polymorphs) could be suppressed.

1 Crystal structure calculation device, 2 Initial structure generation unit, 3 Lattice energy calculation unit, 4, 6, 23 Minimal condition determination unit, 7 Output unit, 8 Atomic position variable change unit, 9 Lattice vector variable change unit, 10, 22 Variables Change unit, 11 control unit, 12 operation unit, 13 storage unit, 21 deviation direction calculation unit, 31 parameter storage unit, 41 empirical expression storage unit, 51 X-ray data storage unit, T1 reciprocal lattice set, T2 subset, V1 to V6 , V1a, V2a, V1'to V3' Lattice vector variable and atomic position variable set, L contour line, U1 to U9 lattice energy, R1 to R9 crystal structure candidates, f free energy continuous function, N1 to N3 deviation direction.

Claims (7)

The initial crystal structure of the target material including the reciprocal lattice vector indicating the reciprocal lattice vector of the unit cell and the atomic position variable indicating the position vector of each atom in the unit cell is obtained by the lattice vector variable of the target material and the atomic position. Change sequentially within the range of the reciprocal lattice with the variable,
The free energy of the modified set of the lattice vector variable and the atomic position variable was calculated, and it was calculated in a continuous function of free energy whose domain is a continuous subset satisfying the following equation among the Cartesian products. A crystal structure calculation method for calculating a set of the lattice vector variable and the atomic position variable of the target substance based on whether or not the free energy shows a minimum.
Here, _Ust is the lattice energy of the crystal, the vectors a, b and c are the lattice vector variables, the vector r _i (i = 1, 2, ..., N) is the atomic position variable, and N is the atom in the unit cell. Indicates the number of. The lattice vector variable is sequentially changed as an independent variable within the range of the direct product set.
By calculating the atomic position variable whose lattice energy is the minimum with respect to the changed lattice vector variable as the dependent variable, the pair of the lattice vector variable and the atomic position variable included in the subset is searched. ,
A claim that sequentially calculates the free energy of a set of the lattice vector variable and the atomic position variable included in the subset to calculate the set of the lattice vector variable and the atomic position variable indicating the minimum free energy. The crystal structure calculation method according to 1. The lattice vector variable and the atomic position variable are sequentially changed as independent variables within the range of the Cartesian product.
Determines whether the modified lattice vector variable and atomic position variable pair are included in the subset and have the smallest free energy in the subset, based on Lagrange's method of undetermined multipliers. The crystal structure calculation method according to claim 1. The crystal structure calculation according to any one of claims 1 to 3, wherein the free energy is calculated based on the natural frequency of phonons calculated by pseudo-harmonic oscillator approximation from quantum chemistry calculation or classical force field calculation. Method. An empirical formula for calculating the free energy is calculated in advance based on the values of the lattice vector variable and the atomic position variable included in the crystal structure information obtained in the past, and the free energy is calculated based on the empirical formula. The crystal structure calculation method according to any one of claims 1 to 3 to be calculated. A crystal structure calculation program for causing a computer to execute each step of the crystal structure calculation method according to any one of claims 1 to 5. The initial crystal structure of the target substance including the reciprocal lattice vector indicating the reciprocal lattice vector of the unit lattice and the atomic position variable indicating the position vector of each atom in the unit lattice is the direct product of the lattice vector variable and the atomic position variable. A variable change part that changes sequentially within the range of the set,
The free energy of the modified set of the lattice vector variable and the atomic position variable was calculated, and it was calculated in a continuous function of free energy whose domain is a continuous subset satisfying the following equation among the Cartesian products. A crystal structure calculation device including a minimum condition determination unit that calculates a set of the lattice vector variable and the atomic position variable of the target substance based on whether or not the free energy shows a minimum.
Here, _Ust is the lattice energy of the crystal, the vectors a, b and c are the lattice vector variables, the vector r _i (i = 1, 2, ..., N) is the atomic position variable, and N is the atom in the unit cell. Indicates the number of.