CP2K Open Source Molecular Dynamics

User Tools

Site Tools

Trace:
exercises:2015_pitt:gga

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
exercises:2015_pitt:gga [2015/02/26 12:56] – vondeleexercises:2015_pitt:gga [2020/08/21 10:15] (current) – external edit 127.0.0.1
Line 23: Line 23:
   * Use vmd to vizualize the geometries (provided below) named ''mode1.xyz'' and ''mode2.xyz''   * Use vmd to vizualize the geometries (provided below) named ''mode1.xyz'' and ''mode2.xyz''
   * To edit the input files provided below, use an editor such as ''vi'' or ''nano''. While ''nano'' is simple to use, ''vi'' can be [[http://www.cp2k.org/tools| configured to colour-code cp2k inputs]].   * To edit the input files provided below, use an editor such as ''vi'' or ''nano''. While ''nano'' is simple to use, ''vi'' can be [[http://www.cp2k.org/tools| configured to colour-code cp2k inputs]].
-  * You will need files named ''BASIS_MOLOPT'' ''GTH_POTENTIALS'' ''dftd3.dat'' that are provided as part of CP2K in a directory ''cp2k/data'' +  * You will need files named ''BASIS_MOLOPT'' ''GTH_POTENTIALS'' ''dftd3.dat'' that are provided as part of CP2K in a directory ''cp2k/data'', unless the code has been compiled with the proper flag ('' -D__DATA_DIR'') so that these are found automatically. 
-  * Use a job script to submit jobs on the cluster, an example job submission script might look like (TODO: adjust for environment)+  * Use a job script to submit jobs on the cluster, an example job submission script might look like 
  
 <code - job> <code - job>
-#!/bin/bash --login +#PBS -N mode1 
-#SBATCH --job-name=mode1 +#PBS -j oe 
-#SBATCH --nodes=6 +#PBS -q dist_small 
-#SBATCH --time=0:10:00 +#PBS -l nodes=4:ppn=16 
-#SBATCH --account=ch5 +#PBS -l walltime=10:00 
-aprun -n 144 -N 24 -d 1 ../../cp2k.popt -i mode1.inp -o mode1.out +#PBS -A cp2k2015
-</code>+
  
 +cd $PBS_O_WORKDIR
  
-===== 2. Task: Bond induced density differences ===== +module purge 
-Compute the density difference induced by the bonding for the first binding mode. +module load cp2k/2.6
-For this you will have to run three separate energy calculations: +
-  - combined system bound in the first mode (file ''mode1.xyz'') +
-  - lone acetic acid molecule (just remove slab's coordinates from ''mode1.xyz'') +
-  - lone TiO$_2$ slab (just remove the acid's coordinates from ''mode1.xyz'')+
  
-In order to output the electronic densities as cube files, your input file has to contain the following snipped: +prun cp2k.popt -i mode1.inp -o mode1.out
-<code> +
-&DFT +
-  &PRINT +
-    &E_DENSITY_CUBE +
-    &END E_DENSITY_CUBE +
-  &END +
-&END DFT+
 </code> </code>
  
-<note tip> 
-The calculations involving the large TiO$_2$ slab should be run on 16 nodes with ''bsub  -n 16''. 
-</note> 
  
-To process the cube files we are going to use the cubecruncher tool. It is part of CP2K, but not installed on brutus. +===== 2. Task: Binding induced density differences =====
-Therefore, a compiled binary of the tool is provided at ''/cluster/home03/matl/schuetto/cubecruncher.x''. Before invoking the cubecruncher, you have to load the cp2k module. To simplify the following steps you should create a symbolic link to it in your working directory: +
-<code> +
-you@brutusX ~$ ln -s /cluster/home03/matl/schuetto/cubecruncher.x . +
-</code>+
  
 +We start with single point energy calculations on binding mode 1, to visualize the interaction between molecule and surface. The goal is to compute the binding induced density difference:
 +
 +\[ \rho_\text{induced}= \rho_\text{slab-dye-complex} - \rho_\text{dye} - \rho_\text{slab} \]
 +
 +First, we'll discuss in detail the structure and the choices made in the sample input file ''mode1.inp''.
 +
 +topics:
 +  * Project name
 +  * Runtype
 +  * Gaussian Basis, pseudopotentials
 +  * PW Cutoff
 +  * thresholds
 +  * SCF: OT
 +  * XC and -D3 correction
 +  * Unit cell choice
 +
 +Second, we run the cp2k input and store the output for analysis and discussion.
  
-Now subtract the densities of the lone systems from the bonded system: 
 <code> <code>
-you@brutusX ~$ ./cubecruncher.x -i mode1-ELECTRON_DENSITY-1_0.cube -subtract mode1_acid-ELECTRON_DENSITY-1_0.cube -o tmp.cube +cp2k.popt -i mode1.inp -o mode1.out
-you@brutusX ~$ ./cubecruncher.x -i tmp.cube -subtract mode1_slab-ELECTRON_DENSITY-1_0.cube -o mode1_delta.cube+
 </code> </code>
  
 +In addition to the output ''mode1.out'', other files are gerenated by CP2K named MODE1*
  
-The generated cube file is not aligned with the simulation cell. Center the cube file with the cubecruncher.x tool:+topics: 
 +  * General overiew 
 +  * OT output 
 +  * Various grid quantities 
 +  * Density cube output 
 +  * Timing report 
 + 
 +Third, we compute the changes in density induced by the binding. For this you will have to run three separate energy calculations: 
 +  - combined system bound in the first mode (file ''mode1.xyz'') 
 +  - lone acetic acid molecule (just remove slab's coordinates from ''mode1.xyz''), name the file ''mode1_dye.xyz'' 
 +  - lone TiO$_2$ slab (just remove the acid's coordinates from ''mode1.xyz''), name the file ''mode1_slab.xyz'' 
 + 
 +Create the .xyz files (check with vmd that they contain the right subsystems), and create mode1_dye.inp, mode1_slab.inp by changing both ''COORD_FILE_NAME'' and ''PROJECT'' accordingly. 
 + 
 +After computing these input files, we analyze the results using a tool provided with cp2k ''cubecruncher.x'' that manipulates cube files (e.g. can compute difference).
 <code> <code>
-you@brutusX ~$ cubecruncher.x -center geo -i mode1_delta.cube -o mode1_delta-centered.cube+~$ cubecruncher.x -i MODE1-ELECTRON_DENSITY-1_0.cube -subtract MODE1_dye-ELECTRON_DENSITY-1_0.cube -o tmp.cube 
 +~$ cubecruncher.x -i tmp.cube -subtract MODE1_slab-ELECTRON_DENSITY-1_0.cube -center geo -o MODE1_delta.cube
 </code> </code>
  
-You can visualize the resulting file ''mode1_delta-centered.cube'' with VMD. This has been covered in a [[mo_ethene| previous exercise]].+You can visualize the resulting file ''MODE1_delta.cube'' with VMD. 
  
 What you get should look similar to this:  What you get should look similar to this: 
-{{ dye_tio_bonding_density.png?300 |}}+{{ exercises:2015_ethz_mmm:dye_tio_bonding_density.png?300 |}} 
 + 
 +===== 3. Task: relative stabilities  ===== 
 + 
 +In order to compute the relative stability of mode1 and mode2, both configurations need to be geometry optimized.  
 + 
 +To do so, turn off the generation of cubes (''&E_DENSITY_CUBE OFF'') in mode1.inp, change to ''RUN_TYPE GEO_OPT'' and adjust the project name. Create and run a similar input file for mode2. 
 + 
 +input topics: 
 +  * BFGS vs LBFGS 
 +  * EPS_SCF, CUTOFF, MAX_DR, .. 
 + 
 +output topics: 
 +  * ''Informations at step'' 
 +  * Trajectory ''MODE1_GEO-pos-1.xyz'' 
 + 
 +Compare the final energies (''ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):''), and determine which mode is most stable. Does this agree with the values in table 1 of the manuscript cited ? 
 +===== 4. Task: ab initio molecular dynamics  =====
  
-===== 3. Task: Bonding energies  ===== +<note>This task is optional, and can be done near the end if time is available. Ab initio MD will be briefly discussed in a next exercise.</note>
-Compute the binding energy for both binding modes:+
  
-\[ E_\text{binding}=\sum E_\text{products} - \sum E_\text{reactants} \]+Perform a short ab initio molecular dynamics simulation of the system (~1000 steps, ~0.5ps) by changing to ''RUN_TYPE MD''. After a couple of hours the job should be finished. Now analyze the OH distance in VMD. A possible outcome is :
  
-For this you will need the energy values of four systems: +{{exercises:2015_pitt:acetic_acid_mode1_md.png?600 |}}
-  - lone acetic acid molecule (run geometry optimization, use energy of last step) +
-  - lone TiO$_2$ slab  (you can use the already geometry optimized coordinates from ''relaxed_slab.xyz'' at the end of the exercise) +
-  - combined system bound in the first mode (can be reused from previous task) +
-  - combined system bound in the second mode (file ''mode2.xyz'')+
  
-<note important> +What can you say about the hydrogen bond to the surface, relative acidity of the two oxygens ?  
-You can not reuse the energy values for the lone sub-systems from the previous task. Since the unbound subsystems might relax into a different geometry, they have to be geometry optimized first. This has been covered in a  +Note that, in order to be statistically relevant, longer trajectories should be employed, and surface slab thickness will play an important role. Also compare to Fig. 7 of the paper referenced. 
-[[geometry_optimization|previous exercise]]. +====== Required Files ======
-</note>+
  
 +(right) click on the filename to download to your local machine.
  
-===== Required Files ===== 
 <code - mode1.inp> <code - mode1.inp>
 &GLOBAL &GLOBAL
Line 121: Line 145:
     &MGRID     &MGRID
        ! PW cutoff ... depends on the element (basis) too small cutoffs lead to the eggbox effect.        ! PW cutoff ... depends on the element (basis) too small cutoffs lead to the eggbox effect.
-       ! certain calculations (e.g. geometry optimization, vibrational frequencies, NPT and cell optimizations, need higher cutoffs)+       ! certain calculations (e.g. geometry optimization, vibrational frequencies, 
 +       ! NPT and cell optimizations, need higher cutoffs)
        CUTOFF [Ry] 400         CUTOFF [Ry] 400 
     &END     &END
  
     &QS     &QS
-       METHOD GPW  ! use the GPW method (i.e. pseudopotential based calculations with the Gaussian and Plane Waves scheme). +       ! use the GPW method (i.e. pseudopotential based calculations with the Gaussian and Plane Waves scheme). 
-       EPS_DEFAULT 1.0E-10 ! default threshold for numerics ~ roughly numerical accuracy of the total energy per electron, sets reasonable values for all other thresholds +       METHOD GPW  
-       EXTRAPOLATION ASPC ! used for MD, the method used to generate the initial guess. +       ! default threshold for numerics ~ roughly numerical accuracy of the total energy per electron, 
 +       ! sets reasonable values for all other thresholds. 
 +       EPS_DEFAULT 1.0E-10  
 +       ! used for MD, the method used to generate the initial guess. 
 +       EXTRAPOLATION ASPC 
     &END     &END
  
Line 185: Line 214:
     &END CELL     &END CELL
  
-    ! atom coordinates can be in the &COORD section, or provided as an external file.+    ! atom coordinates can be in the &COORD section, 
 +    ! or provided as an external file.
     &TOPOLOGY     &TOPOLOGY
       COORD_FILE_NAME mode1.xyz       COORD_FILE_NAME mode1.xyz
Line 191: Line 221:
     &END     &END
  
-    ! MOLOPT basis sets are fairly costly, but in the 'DZVP-MOLOPT-SR-GTH' available for all elements +    ! MOLOPT basis sets are fairly costly, 
-    ! their contracted nature makes them suitable for condensed and gas phase systems alike.+    ! but in the 'DZVP-MOLOPT-SR-GTH' available for all elements 
 +    ! their contracted nature makes them suitable 
 +    ! for condensed and gas phase systems alike.
     &KIND H                                   &KIND H                              
       BASIS_SET DZVP-MOLOPT-SR-GTH                BASIS_SET DZVP-MOLOPT-SR-GTH         
Line 225: Line 257:
    TEMPERATURE [K] 300    TEMPERATURE [K] 300
    TIMESTEP [fs] 0.5    TIMESTEP [fs] 0.5
-   STEPS 100+   STEPS 1000
    # GLE thermostat as generated at http://epfl-cosmo.github.io/gle4md     # GLE thermostat as generated at http://epfl-cosmo.github.io/gle4md 
    # GLE provides an effective NVT sampling.    # GLE provides an effective NVT sampling.
exercises/2015_pitt/gga.1424955413.txt.gz · Last modified: 2020/08/21 10:14 (external edit)