:sequential_nav: next .. _tutorial-basic-vibrational-properties-benzene: Modes of vibration of the benzene molecule ========================================== .. :Author: Javier Junquera and Andrei Postnikov .. sidebar:: **Have you set up the local environment?** If not, :ref:`do that now ` before proceeding. The goal of this exercise is to compute the vibrational frequency of a molecule (benzene). .. hint:: Enter the ``Benzene`` directory. Before running the calculation to compute the vibrational frequencies, the first step is to relax the geometrical structure of the system under study. To start with, we run a conjugate gradient minimization to relax the atomic positions. The input file has been prepared for you in the file benzene.relax.fdf. See how the structure of the benzene molecule has been introduced in the **Z-matrix format** (especification of internal variables, such as distances, angles, and torsional angles). This allows the minimization including some constraints in the symmetry in a trivial way. We just allow the C-C and C-H distances to relax:: %block Zmatrix molecule 2 0 0 0 xm1 ym1 zm1 0 0 0 2 1 0 0 CC 90.0 60.0 0 0 0 2 2 1 0 CC CCC 90.0 0 0 0 2 3 2 1 CC CCC 0.0 0 0 0 2 4 3 2 CC CCC 0.0 0 0 0 2 5 4 3 CC CCC 0.0 0 0 0 1 1 2 3 CH CCH 180.0 0 0 0 1 2 1 7 CH CCH 0.0 0 0 0 1 3 2 8 CH CCH 0.0 0 0 0 1 4 3 9 CH CCH 0.0 0 0 0 1 5 4 10 CH CCH 0.0 0 0 0 1 6 5 11 CH CCH 0.0 0 0 0 constants ym1 5.00 zm1 0.00 CCC 120.0 CCH 120.0 variables CC 1.390 CH 1.090 constraints xm1 CC -1.0 3.903229 %endblock Zmatrix We run Siesta to carry out the relaxation:: siesta benzene.relax.fdf > benzene.relax.out Computing the force constants in real space -------------------------------------------- We have prepared an input file, called benzene.ifc.fdf to run Siesta and compute the interatomic force constants in real space. In this, we have already copied the relaxed coordinates and unit cell from the *benzene.XV* generated previously into the `AtomicCoordinatesAndAtomicSpecies` block. We have also included the atomic masses after the coordinates of each atom. This will be useful for the next step with the *vibra* code Here are the relevant sections:: LatticeConstant 1.0 Bohr %block LatticeVectors 21.938124322 0.000000000 0.000000000 0.000000000 20.556799916 0.000000000 0.000000000 0.000000000 11.910755412 %endblock LatticeVectors AtomicCoordinatesFormat NotScaledCartesianBohr %block AtomicCoordinatesAndAtomicSpecies 4.732644349 9.448630623 0.000000000 2 12.01070 6.054339080 11.737873050 0.000000000 2 12.01070 8.697728543 11.737873050 0.000000000 2 12.01070 10.019423274 9.448630623 0.000000000 2 12.01070 8.697728543 7.159388196 0.000000000 2 12.01070 6.054339080 7.159388196 0.000000000 2 12.01070 2.637688094 9.448630623 0.000000000 1 1.00794 5.006860953 13.552158387 0.000000000 1 1.00794 9.745206671 13.552158387 0.000000000 1 1.00794 12.114379529 9.448630623 0.000000000 1 1.00794 9.745206671 5.345102860 0.000000000 1 1.00794 5.006860953 5.345102806 0.000000000 1 1.00794 %endblock AtomicCoordinatesAndAtomicSpecies .. note:: If you want to copy your own coordinates from the XV file, you'll need to do some rearrangements. Remember that for vibra, the ``AtomicCoordinatesAndAtomicSpecies`` block needs to look like this:: X Y Z SpeciesIndex Mass To compute the interatomic force constant in real space, we have to run Siesta:: siesta < benzene.ifc.fdf > benzene.ifc.out The interatomic force constant matrix in real space is stored in a file with name of the form *SystemLabel.FC*. .. Again, the explanation of the different entries of this file can be found in the theoretical lectures and in the manual. Computing the dynamical matrix at the Gamma point, and (phonon) modes --------------------------------------------------------------------- Once the interatomic force constants in real space have been computed, a discrete Fourier transform is performed to compute the dynamical matrix in reciprocal space. Then, the dynamical matrix is diagonalized and its eigenvalues and eigenvectors are computed. This is done using the VIBRA code. In the case of a molecule, only the Gamma point is relevant. It is specified in the same way as to compute the electronic band structure, in the same file benzene.ifc.fdf:: Eigenvectors .true. # Compute both phonon eigenvalues and eigenvectors BandLinesScale pi/a %block BandLines 1 0.0 0.0 0.0 \Gamma # Only the Gamma point (enough for a molecule) %endblock BandLines To compute the vibrational frequencies:: vibra < benzene.ifc.fdf > vibra.out The output of this code is: * *SystemLabel.bands*: with the different mode frequencies (in cm^-1). They are stored in the same way as the electronic band structure. * *SystemLabel.vectors*: with the eigenmodes at Gamma (the format is self-explained). .. note:: Some of the modes might have negative frequencies. How could that be? How to visualize the normal modes --------------------------------- After getting the .vectors file (calculated by vibra) and the .XV file (computed in Siesta), run the vib2xsf program. You have to answer a few question on the fly, regarding the name of the files where the .vectors are stored, the units to be used to introduce the lattice vectors (Bohrs or Angstroms), the zero of coordinates, the unit cell lattice vectors, the first mode to visualize, the last mode to visualize, the amplitude of the modes to be visualized, and the number of steps in the movie. You can play a little bit, but to save time we have prepared all the answers in the file vib2xsf.dat for you. Just run:: vib2xsf < vib2xsf.dat This will produce two files per mode: * .XSF file: contains a static structures (as in .XV) * .AXSF file: contains the animation of a phonon, for a (user-chosen) amplitude and number of steps. They can be visualized using OVITO:: $ ovito benzene.Mode_1.AXSF (Animation XCrySDen Structure File) Or open ovito :: $ ovito File > Load File Load benzen.Mode_*.AXSF .. note:: It might be interesting to analyze those modes with negative or small frequencies.