This is an old revision of the document!
Table of Contents
Run your first simulation using CP2K
When you check CP2K's features and the outline of the lecture you will notice that there are many levels of theory, methods and possibilities to combine them. This results in a large number of possible options and coefficients to setup, control and tune a specific simulation. Together with the parameters for the numerical solvers this means that an average CP2K configuration file will contain quiet a number of options (even though for many others the default value will be applied) and not all of them will be discussed in the lecture or the exercises.
The CP2K Manual is the complete reference for all configuration options. Where appropriate you will find a reference to the respective paper when looking up a specific keyword/option.
To get you started, we will do a simple exercise using Molecular Mechanics (that is: a classical approach). The point is to get familiar with the options, organizing and editing the input file and analyze the output.
Computation of the Lennard Jones curve
In this exercise you will compute the Lennard-Jones energy curve for a system of two Krypton (Kr) atoms using a Molecular Mechanics simulation rather than the analytical form of the potential.
In Part I you find the instructions for computing the energy of two Kr atoms at a distance $r=4.00 Å$.
In Part II you find the instructions for getting the energy profile as a function of $r$.
Additonal parameters for Neon (Ne) and combination rules to obtain new parameters are provided in Part III and IV.
You are expected to hand in the respective plots plus answers to the questions. The format can be either a proper PDF or directly in the email with the plots attached (as pictures in png, jpg, or pdf format).
Part I: Single Point (Energy) calculation
In this section a commented CP2K input example for a single point calculation is provided. Comments are added, they start with an exclamation mark '!'.
1. Step
Load the CP2K module as shown before, create a directory ex1 and change to it:
mkdir ex1 cd ex1
Save the following input to a file named energy.inp (for example using $ vim energy.inp):
- energy.inp
&GLOBAL ! section to select the kind of calculation RUN_TYPE ENERGY ! select type of calculation. In this case: ENERGY (=Single point calculation) &END GLOBAL &FORCE_EVAL ! section with parameters and system description METHOD FIST ! Molecular Mechanics method &MM ! specification of MM parameters &FORCEFIELD ! parameters needed to describe the potential &SPLINE EMAX_SPLINE 10000 ! numeric parameter to ensure calculation stability. Should not be changed &END &NONBONDED ! parameters for the non bonded interactions &LENNARD-JONES ! Lennard-Jones parameters atoms Kr Kr EPSILON [K_e] 164.56 SIGMA [angstrom] 3.601 RCUT [angstrom] 25.0 &END LENNARD-JONES &END NONBONDED &CHARGE ATOM Kr CHARGE 0.0 &END CHARGE &END FORCEFIELD &POISSON ! solver for non periodic calculations PERIODIC NONE &EWALD EWALD_TYPE none &END EWALD &END POISSON &END MM &SUBSYS ! system description &CELL ABC [angstrom] 10 10 10 PERIODIC NONE &END CELL &COORD UNIT angstrom Kr 0 0 0 Kr 4 0 0 &END COORD &END SUBSYS &END FORCE_EVAL
2. Step
Run a CP2K calculation with the following command:
$ cp2k.sopt -i energy.inp -o energy.out
$ cp2k.sopt -i energy.inp | tee energy.out
which will write the output simultaneously to a file energy.out and show it in the terminal.
3. Step
The output ($ less energy.out) should look like this:
[...] **** **** ****** ** PROGRAM STARTED AT 2016-09-22 15:15:15.977 ***** ** *** *** ** PROGRAM STARTED ON tcopt3 ** **** ****** PROGRAM STARTED BY studentXX ***** ** ** ** ** PROGRAM PROCESS ID 112277 **** ** ******* ** PROGRAM STARTED IN /data/students/studentXX/ex1 [...] ENERGY| Total FORCE_EVAL ( FIST ) energy (a.u.): 0.003617048870059 [...] The number of warnings for this run is : 0 ------------------------------------------------------------------------------- **** **** ****** ** PROGRAM ENDED AT 2016-09-22 15:15:16.027 ***** ** *** *** ** PROGRAM RAN ON tcopt3 ** **** ****** PROGRAM RAN BY studentXX ***** ** ** ** ** PROGRAM PROCESS ID 112277 **** ** ******* ** PROGRAM STOPPED IN /data/students/studentXX/ex1
If you get the closing banner you know that CP2K finished.
The number of warnings for this run is : ...
line. If that number is not zero you must check the rest of the output for warnings and act accordingly.
The following line tells you the result:
ENERGY| Total FORCE_EVAL ( FIST ) energy (a.u.): -0.000518941408898
This is the energy (in Hartree) for a system of 2 Kr atoms at distance $ r=4.00 Å$
Note, that in the input-file EPSILON is given in units of Kelvin, whereas in the output the energy is printed in Hartree, which is the unit of energy in the system of atomic units (a.u.).
To convert from Kelvin to Hartree you have to multiply with the Boltzmann constant $ k_\text{b} = 3.1668154 \cdot 10^{-6} \frac{E_\text{H}}{\text{K}} $ .
Part II: Computation of the LJ energy curve
In this section a few scripts to get the LJ energy profiles are presented.
1. Step
In order to get a good profile, a set of energy values as a function of the interatomic distance is needed. You can use the energy.inp input file and change the Kr coordinates in order to get different starting distances.
To do so:
$ mv energy.out energy_dist4A.out
For doing so:
$ cp energy.inp energy_dist2A.inp
then edit the input file to update to the new coordinates (e.g. 2 Å) and rerun CP2K to produce a new output file:
$ cp2k.sopt -i energy_dist2A.inp -o energy_dist2A.out
2. Step
When you have tested a few distances, you can produce a table looking like:
| Input file | Distance (Å) | Energy (Eh) |
|---|---|---|
| energy_dist1A.inp | 1 | … |
| energy_dist1.5A.inp | 1.5 | … |
| energy_dist2A.inp | 2 | … |
| energy_dist2.5A.inp | 2.5 | … |
| energy_dist3A.inp | 3 | … |
| … | … | … |
This is the Lennard Jones energy curve for two Kr atoms.
By using any plotting program you can now get a representation of the energy profile.
3. Step
Here are reported the LJ parameters for Ne atoms. Those are to replace the Kr parameters in the input file, along with your Ne coordinates that have to replace the Kr coordinates. A new LJ curve for Ne atoms can be now generated.
&NONBONDED
&LENNARD-JONES ! Lennard-Jones Ne parameters
atoms Ne Ne
EPSILON [K_e] 36.831
SIGMA [angstrom] 2.775
RCUT [angstrom] 25.0
&END LENNARD-JONES
&END NONBONDED
&CHARGE
ATOM Ne
CHARGE 0.0
&END CHARGE
4. Step
Here are reported the combination rules for pairs unlike pairs, i.e. for pairs of non identical atoms.
Once generated the ε and σ parameters for the couple Kr/Ne, generate once more the LJ dissociation curve.
Compare the “mixed” curve to the two “pure” curves and report the position and depth of the minimum.
$$ \sigma_{ij}= \sqrt{\sigma_i\sigma_j}$$
$$ \epsilon_{ij}= \sqrt{\epsilon_i\epsilon_j}$$
- The “ LENNARD-JONES ” section must be present for all the three possible couples: Kr-Kr, Ne-Ne and Ne-Kr. Example:
&LENNARD-JONES ! Lennard-Jones parameters for Ar-Ar interaction
atoms Kr Kr
EPSILON [K_e] 164.56
SIGMA [angstrom] 3.601
RCUT [angstrom] 25.0
&END LENNARD-JONES
&LENNARD-JONES ! Lennard-Jones Ne-Ne parameters
atoms Ne Ne
EPSILON [K_e] 36.831
SIGMA [angstrom] 2.775
RCUT [angstrom] 25.0
&END LENNARD-JONES
&LENNARD-JONES ! Lennard-Jones parameters for Kr-Ne interaction
atoms Kr Ne
EPSILON [K_e] YOUR EPSILON
SIGMA [angstrom] YOUR SIGMA
RCUT [angstrom] 25.0
&END LENNARD-JONES
- The “ CHARGE ” section must be also duplicated:
&CHARGE
ATOM Ne
CHARGE 0.0
&END CHARGE
&CHARGE
ATOM Kr
CHARGE 0.0
&END CHARGE
Questions
- Sketch the LJ energy curve for the two set of parameters ($\sigma$ and $\epsilon$) provided.
- Report, for both curves, the minimum energy distance and the depth of the minimum.
- What are the major differences between the curves? How do they relate to the sets of parameters provided?