exercises:2015_cecam_tutorial:mtd1
Differences
This shows you the differences between two versions of the page.
Both sides previous revisionPrevious revisionNext revision | Previous revision | ||
exercise:simple_metadynamics_simulation_using_the_coordination_numbers_as_variables [2015/07/20 09:01] – correct link to DIRECTION WALL_MINUS tmueller | exercises:2015_cecam_tutorial:mtd1 [2020/08/21 10:15] (current) – external edit 127.0.0.1 | ||
---|---|---|---|
Line 1: | Line 1: | ||
+ | ====== Simple metadynamics simulation using the coordination numbers as variables ====== | ||
+ | |||
Problem: Dissociation reaction of nitric acid on graphene and atomic rearrangements of a < | Problem: Dissociation reaction of nitric acid on graphene and atomic rearrangements of a < | ||
+ | |||
+ | * Original author: Marcella Iannuzzi | ||
+ | * Complete source and output files: [[http:// | ||
+ | |||
+ | ===== Introduction ===== | ||
For this tutorial some input and output files are given in order to present a complete procedure to solve the given problem. Some hints are also given to help in the analysis of the results. In order to be able to run these examples, some paths need to be correctly set in the input files (i.e. set the variables '' | For this tutorial some input and output files are given in order to present a complete procedure to solve the given problem. Some hints are also given to help in the analysis of the results. In order to be able to run these examples, some paths need to be correctly set in the input files (i.e. set the variables '' | ||
Line 8: | Line 15: | ||
* Set up and preliminary simulations to learn about the dynamics of nitric acid on graphene, as obtained at the DFTB level of theory | * Set up and preliminary simulations to learn about the dynamics of nitric acid on graphene, as obtained at the DFTB level of theory | ||
- | * Mmetadynamics | + | * Metadynamics |
* Set up and preliminary simulations to learn about the dynamics of the selected small Si cluster saturated by H atoms | * Set up and preliminary simulations to learn about the dynamics of the selected small Si cluster saturated by H atoms | ||
- | * Metadynamics simulation aimed at observing atomic rearrangements of the cluster bay changing the coordination of both Si and H species. | + | * Metadynamics simulation aimed at observing atomic rearrangements of the cluster bay changing the coordination of both Si and H species. |
| | ||
- | ====== First task: dynamics of two HNO3 molecules over a graphene sheet ====== | + | ===== First task: dynamics of two HNO3 molecules over a graphene sheet ===== |
| | ||
The examples on this system are in the directory '' | The examples on this system are in the directory '' | ||
However, the reaction should be catalyzed in presence of C particles (sooth). In the proposed example, the molecules are located in the vicinity of graphene, that should mimic the role of sooth. Graphene indeed is not very reactive, better models can be considered using defective or functionalized graphene. | However, the reaction should be catalyzed in presence of C particles (sooth). In the proposed example, the molecules are located in the vicinity of graphene, that should mimic the role of sooth. Graphene indeed is not very reactive, better models can be considered using defective or functionalized graphene. | ||
| | ||
- | The initial coordinates are given in '' | + | The initial coordinates are given in '' |
< | < | ||
Line 54: | Line 61: | ||
&END | &END | ||
</ | </ | ||
- | Otherwise the '' | + | Otherwise the '' |
- | The output of this short simulation (5 ps) is stored in '' | + | The output of this short simulation (5 ps) is stored in '' |
< | < | ||
Line 78: | Line 85: | ||
The next step is to set up a few collective variables (CV) that can be later used for metadynamics (MTD) simulations. It is important to select good CV that can describe the relevant configuration along the reaction path. Moreover, it is useful to learn about the typical behavior of the selected CV along an unbiased MD run. Hence, after selecting a set of CV, preliminary runs should be performed in order to monitor the dynamics of these variables. This can be done setting up MTD simulations, | The next step is to set up a few collective variables (CV) that can be later used for metadynamics (MTD) simulations. It is important to select good CV that can describe the relevant configuration along the reaction path. Moreover, it is useful to learn about the typical behavior of the selected CV along an unbiased MD run. Hence, after selecting a set of CV, preliminary runs should be performed in order to monitor the dynamics of these variables. This can be done setting up MTD simulations, | ||
- | The input '' | + | The input '' |
< | < | ||
Line 147: | Line 154: | ||
</ | </ | ||
- | This last CV controls the distance of the molecules from the layer, which is an important factor to determine whether the dissociation is somehow favored by the presence of graphene. The output of the unbiased MD that monitors the behavior of these 4 CV is in '' | + | This last CV controls the distance of the molecules from the layer, which is an important factor to determine whether the dissociation is somehow favored by the presence of graphene. The output of the unbiased MD that monitors the behavior of these 4 CV is in '' |
< | < | ||
Line 186: | Line 193: | ||
| | ||
</ | </ | ||
- | With '' | + | With '' |
By plotting the CV as recorded along the short MD trajectory (3 ps), the amplitude of the equilibrium fluctuations can be evaluated and then used to set up the size of the Gaussian hills that build up the biasing potential. The first CV fluctuates close to zero, with fluctuations smaller than 0.2. The second is around 2.8. The fluctuations are smaller due to the stiffness of the three NO bonds. The coordination of H to C is also typically zero, but it can change a lot when the molecules approach the layer, even if there is no dissociation of H and no binding to C. This indicates that this variable is difficult to control and might turn out to be tricky to use it to distinguish among different states of the reaction process. The point to plane distance shows quite large fluctuations and it is clearly not suited to distinguish a specific state along the reaction path. Moreover, its minima, when the two molecules are closer to the layer, correspond to the maxima of the third CV, i.e. the CN of H to C. At least before dissociation, | By plotting the CV as recorded along the short MD trajectory (3 ps), the amplitude of the equilibrium fluctuations can be evaluated and then used to set up the size of the Gaussian hills that build up the biasing potential. The first CV fluctuates close to zero, with fluctuations smaller than 0.2. The second is around 2.8. The fluctuations are smaller due to the stiffness of the three NO bonds. The coordination of H to C is also typically zero, but it can change a lot when the molecules approach the layer, even if there is no dissociation of H and no binding to C. This indicates that this variable is difficult to control and might turn out to be tricky to use it to distinguish among different states of the reaction process. The point to plane distance shows quite large fluctuations and it is clearly not suited to distinguish a specific state along the reaction path. Moreover, its minima, when the two molecules are closer to the layer, correspond to the maxima of the third CV, i.e. the CN of H to C. At least before dissociation, | ||
- | ====== Second task: Metadynamics of the dissociation of HNO3 over a graphene sheet ====== | + | ===== Second task: Metadynamics of the dissociation of HNO3 over a graphene sheet ===== |
- | The presented MTD run employs as CV only the three CN described above. The related input file is '' | + | The presented MTD run employs as CV only the three CN described above. The related input file is '' |
< | < | ||
Line 233: | Line 240: | ||
One Gaussian hill is spawned every NT_HILLS MD steps, while the height of the hill in hartree is given by WW. These parameters together with the width of the Gaussian hills are important to determine the accuracy of the description of the FES through the MTD biasing potential. Since each variable has in principle different dimensions and different dynamics, the shape of the hills filling up the $N_{\text{CV}}$-dimensional configurations space, as defined by the selected CVs, is not isotropic. The parameter SCALE associated to the $i$-th MTD variable determines the amplitude of the Gaussian in the $i$-th space-direction of the $N_{\text{CV}}$-dimensional configuration space. This parameter, as well as the hill’s height and the frequency of collocation, | One Gaussian hill is spawned every NT_HILLS MD steps, while the height of the hill in hartree is given by WW. These parameters together with the width of the Gaussian hills are important to determine the accuracy of the description of the FES through the MTD biasing potential. Since each variable has in principle different dimensions and different dynamics, the shape of the hills filling up the $N_{\text{CV}}$-dimensional configurations space, as defined by the selected CVs, is not isotropic. The parameter SCALE associated to the $i$-th MTD variable determines the amplitude of the Gaussian in the $i$-th space-direction of the $N_{\text{CV}}$-dimensional configuration space. This parameter, as well as the hill’s height and the frequency of collocation, | ||
- | The '' | + | The '' |
The provided trajectory is about 100 ps long and indeed it shows the dissociation of the two molecules into < | The provided trajectory is about 100 ps long and indeed it shows the dissociation of the two molecules into < | ||
- | Other quantities that can be monitored from the '' | + | Other quantities that can be monitored from the '' |
- | ====== Third task: dynamics of Si6H8 ====== | + | ===== Third task: dynamics of Si6H8 ===== |
The data file for this example are in '' | The data file for this example are in '' | ||
- | As in the previous example, preliminary MD runs are carried out to study the dynamics of the system. The electronic structure is computed at the DFT level, using the PBE functional. A standard constant temperature MD is simulated running the input '' | + | As in the previous example, preliminary MD runs are carried out to study the dynamics of the system. The electronic structure is computed at the DFT level, using the PBE functional. A standard constant temperature MD is simulated running the input '' |
- | Possible changes in stoichiometry (loosing H) and structure are going to induce variations in the coordination shell of the Si atoms. Therefore, three CN are selected as CV : Si to Si, Si to H, and H to H. These variables are monitored over the equilibrium trajectory by running the input '' | + | Possible changes in stoichiometry (loosing H) and structure are going to induce variations in the coordination shell of the Si atoms. Therefore, three CN are selected as CV : Si to Si, Si to H, and H to H. These variables are monitored over the equilibrium trajectory by running the input '' |
< | < | ||
Line 260: | Line 267: | ||
the Si-Si CN oscillates around 2. Actually, in the chair configuration 4 Si are three fold coordinated with neighboring Si atoms and 2 are two fold coordinated. The bond length is about 2.4 Å. By changing the curvature of the function, the average value of the CN can be easily moved towards 3 and it can become more sensitive to fluctuations of the Si-Si bond length. The fluctuations of the two other CN are even smaller. Si-H fluctuates around 1.4 and H-H is very close to zero, since in the initial configuration the H atoms do not see each other. | the Si-Si CN oscillates around 2. Actually, in the chair configuration 4 Si are three fold coordinated with neighboring Si atoms and 2 are two fold coordinated. The bond length is about 2.4 Å. By changing the curvature of the function, the average value of the CN can be easily moved towards 3 and it can become more sensitive to fluctuations of the Si-Si bond length. The fluctuations of the two other CN are even smaller. Si-H fluctuates around 1.4 and H-H is very close to zero, since in the initial configuration the H atoms do not see each other. | ||
- | The output of a second simulation that monitors the same three CVs is in MTD_L_H0, and the corresponding input is '' | + | The output of a second simulation that monitors the same three CVs is in MTD_L_H0, and the corresponding input is '' |
< | < | ||
Line 273: | Line 280: | ||
A temperature is associated to the auxiliary variables and can be controlled by temperature rescaling. The use of thermostats for such few degrees of freedom is questionable. The Lagrangian MTD formalism is used in order to better control the kinetics of the CV. This control is obtained through the coupled to the auxiliary variables, whose dynamics depends on the mass and the temperature, | A temperature is associated to the auxiliary variables and can be controlled by temperature rescaling. The use of thermostats for such few degrees of freedom is questionable. The Lagrangian MTD formalism is used in order to better control the kinetics of the CV. This control is obtained through the coupled to the auxiliary variables, whose dynamics depends on the mass and the temperature, | ||
- | Whenever a '' | + | Whenever a '' |
< | < | ||
Line 286: | Line 293: | ||
</ | </ | ||
- | In this case the potential function is $f(s)=K (s-s_0) ^4$, and it is activated whenever the variable becomes lower ('' | + | In this case the potential function is $f(s)=K (s-s_0) ^4$, and it is activated whenever the variable becomes lower ('' |
When such Lagrangian scheme is used more columns appear in the COLVAR, which contain all the relevant information. The 1st column is always the time in fs, the next $N_{\text{CV}}$ columns are the instantaneous values of the auxiliary variables, followed by $N_{\text{CV}}$ columns where the instantaneous values of the CVs are reported. The next columns report the values of the potential gradients: $N_{\text{CV}}$ columns for the gradients of the harmonic potential, $N_{\text{CV}}$ for the gradients of the MTD biasing potential, and $N_{\text{CV}}$ for the gradient of the WALL potential (these are zeros when the corresponding potential is not activated). The following $N_{\text{CV}}$ are the velocities of the auxiliary variables. Then there are the instantaneous values of the harmonic potential, of the MTD potential, and of the WALL potential. The last column is the temperature of the auxiliary variables. | When such Lagrangian scheme is used more columns appear in the COLVAR, which contain all the relevant information. The 1st column is always the time in fs, the next $N_{\text{CV}}$ columns are the instantaneous values of the auxiliary variables, followed by $N_{\text{CV}}$ columns where the instantaneous values of the CVs are reported. The next columns report the values of the potential gradients: $N_{\text{CV}}$ columns for the gradients of the harmonic potential, $N_{\text{CV}}$ for the gradients of the MTD biasing potential, and $N_{\text{CV}}$ for the gradient of the WALL potential (these are zeros when the corresponding potential is not activated). The following $N_{\text{CV}}$ are the velocities of the auxiliary variables. Then there are the instantaneous values of the harmonic potential, of the MTD potential, and of the WALL potential. The last column is the temperature of the auxiliary variables. | ||
- | The dynamics of the CVs along the two simulations, | + | The dynamics of the CVs along the two simulations, |
- | ====== Fourth task: Lagrangian MTD of the atomic rearrangement of Si6H8 ====== | + | ===== Fourth task: Lagrangian MTD of the atomic rearrangement of Si6H8 ===== |
- | '' | + | '' |
exercises/2015_cecam_tutorial/mtd1.1437382908.txt.gz · Last modified: 2020/08/21 10:14 (external edit)