exercises:2019_uzh_acpc2:ex03
Differences
This shows you the differences between two versions of the page.
Both sides previous revisionPrevious revisionNext revision | Previous revision | ||
exercises:2019_uzh_acpc2:ex03 [2019/05/13 22:38] – keimre | exercises:2019_uzh_acpc2:ex03 [2020/08/21 10:15] (current) – external edit 127.0.0.1 | ||
---|---|---|---|
Line 5: | Line 5: | ||
$$ \text{Cl}^- + \text{CH}_3\text{Cl} \longleftrightarrow \text{Cl}\text{CH}_3 + \text{Cl}^-. $$ | $$ \text{Cl}^- + \text{CH}_3\text{Cl} \longleftrightarrow \text{Cl}\text{CH}_3 + \text{Cl}^-. $$ | ||
- | For energy and force evaluation, we will use the Parameterization Method 6 (PM6), which is a relatively inexpensive semiempirical model for electronic structure. For accurate reaction characterizations, | + | For energy and force evaluation, we will use the Parameterization Method 6 (PM6), which is a relatively inexpensive semiempirical model for electronic structure |
===== NEB activation barrier ===== | ===== NEB activation barrier ===== | ||
Line 13: | Line 13: | ||
=== Geometry optimizations === | === Geometry optimizations === | ||
- | To start the NEB calculation, | + | To start the NEB calculation, |
<code - geo.inp> | <code - geo.inp> | ||
Line 135: | Line 135: | ||
</ | </ | ||
- | In the main output of the NEB calculation, | + | In the main output of the NEB calculation, |
< | < | ||
******************************************************************************* | ******************************************************************************* | ||
Line 151: | Line 151: | ||
| | ||
</ | </ | ||
- | The sections show for every replica geometry along the NEB trajectory, the distance to its neighbors and its energy. The final section corresponds to the converged NEB trajectory. | + | These sections show for every replica geometry along the NEB trajectory, the distance to its neighbors and its energy. The final section corresponds to the converged NEB trajectory. |
< | < | ||
Line 163: | Line 163: | ||
Sampling the free energy surface (FES) of a chemical system is a convenient method to explore various stable conformations and possible reaction pathways. To calculate the FES for complicated systems, advanced sampling methods (such as umbrella sampling, metadynamics, | Sampling the free energy surface (FES) of a chemical system is a convenient method to explore various stable conformations and possible reaction pathways. To calculate the FES for complicated systems, advanced sampling methods (such as umbrella sampling, metadynamics, | ||
- | The FES is a projection of the high-dimensional free energy landscape | + | The FES is a projection of the high-dimensional free energy landscape into a small number, usually |
- | To help the system | + | To help the calculation |
The following CP2K input script runs our MD calculation and prints out the CV values for every step: | The following CP2K input script runs our MD calculation and prints out the CV values for every step: | ||
Line 270: | Line 270: | ||
where $s$ is the set CVs and $P(s)$ is the probability that the system has the set of CV values $s$. | where $s$ is the set CVs and $P(s)$ is the probability that the system has the set of CV values $s$. | ||
- | The following Python script can be used to calculate the FES from the '' | + | The following Python script can be used to calculate the FES from the '' |
< | < | ||
- | TBA | + | import numpy as np |
+ | import matplotlib.pyplot as plt | ||
+ | |||
+ | bohr_2_angstrom = 0.529177 | ||
+ | kb = 8.6173303e-5 # eV * K^-1 | ||
+ | |||
+ | temperature = 1000.0 | ||
+ | colvar_path = " | ||
+ | |||
+ | # Load the colvar file | ||
+ | colvar_raw = np.loadtxt(colvar_path) | ||
+ | |||
+ | # Extract the two CVs | ||
+ | d1 = colvar_raw[:, | ||
+ | d2 = colvar_raw[:, | ||
+ | |||
+ | # Create a 2d histogram corresponding to the CV occurances | ||
+ | cv_hist = np.histogram2d(d1, | ||
+ | |||
+ | # probability from the histogram | ||
+ | prob = cv_hist[0]/ | ||
+ | |||
+ | # Free energy surface | ||
+ | fes = -kb * temperature * np.log(prob) | ||
+ | |||
+ | # Save the image | ||
+ | extent = (np.min(cv_hist[1]), | ||
+ | plt.figure(figsize=(8, | ||
+ | plt.imshow(fes.T, | ||
+ | cbar = plt.colorbar() | ||
+ | cbar.set_label(" | ||
+ | plt.xlabel(" | ||
+ | plt.ylabel(" | ||
+ | plt.savefig(" | ||
+ | plt.close() | ||
</ | </ | ||
+ | |||
+ | Here is an example output for temperature 1000K. We clearly see the two local minima corresponding to the one of the chlorine atoms being covalently bonded (distance 1.8 Å) while the other one is around distance 2.5 Å. | ||
+ | |||
+ | {{ : | ||
< | < | ||
- | * Run the MD calculation for 800K, 1000K, 1200K and 1400K. | + | * Run the MD calculation for 400K, 800K, 1200K and 1600K. (The calculations can a take a while.) |
* Create the corresponding FES plots and discuss the temperature dependence. | * Create the corresponding FES plots and discuss the temperature dependence. | ||
- | * In general, how does potential energy differ from free energy? | + | * In general, how does potential energy differ from free energy? |
</ | </ | ||
- |
exercises/2019_uzh_acpc2/ex03.1557787112.txt.gz · Last modified: 2020/08/21 10:15 (external edit)