Category Archives: Molecular Modelling

IR and Raman spectra calculated with Quantum Espresso

Infrared (IR) and Raman spectra are two very valuable tools for the characterization of chemical compounds. And although there seem to be many different possibilities to produce them computationally, I didn’t find any clear tutorial on how to produce them using free software tools. Some of them, give complicated instructions on how to calculate second derivatives or Raman tensors, or stuff like that, but in the end one wonders: Where are the numbers I have to plot!

…Thereby, here is my post about computing IR and Raman spectra with Quantum Espresso.

For the impatient: most of this post is resumed in my Quantum Espresso example for the PHonon code, which will run for some 15-17 minutes and pop-up the IR/Raman spectra of CO2 and ZnO (Wurtzite).

Continue reading →

cube2xyz: convert cube data to xyz and slize segments/planes

In this post, I introduce my new cube2xyz script for converting cube data into “x,y,z,value” format and optionally creating nice 2D or 3D plots displaying the value of the cube property along a segment or a surface.

.cube files are originally from the Gaussian molecular modeling code, but nowadays can be produced and read by most of the computational chemistry programs. They consist of a lot of blocks of data describing the values of a measured property (charge, spin polarization, electronic density, …) on each point of the space inside the simulation cell.

Continue reading →

shrink_traj v0.4: xyz with unit cell info

One of the limitations of xyz coordinate format is generally the lack of extra information, such as the unit cell dimensions or lattice vectors. This can become annoying when we are working with solids, or we want to check whether our molecule has enough space around, in order avoid interaction with its own periodic replicas.
The new shrink_traj v0.4, includes the possibility to include the Jmol readable unit cell info the second line of the xyz file. This way, we can shrink a trajectory (or keep it as it is), and generate a new one which allows us to visualize the unit cell.

It also contains an option to insert Jmol commands directly into the xyz file, such as “background white”, which enables the possibility to customize the view of the structure/trajectory. Continue reading →

Checking the optimum cutoff for QE

Rise your hand if you have ever gotten a pseudo-potential from somewhere and you started to use it just by checking the recommended values of the cutoffs inside the file.

Are we sure these values are good?
If you are an skeptic scientist and you do not thrust what you don’t do by yourself, now there is a really simple and cheap way for you to test for the optimum cutoff of any pseudopotential with Quantum Espresso: my pw_cutoff checking script.

The ecutwfc is described in the manual as the kinetic energy cutoff for wave functions in Ry. Continue reading →

pw2cellvec: monitor vc-relax convergence in QE

When running a cell relaxation, it is sometimes nice to monitor what is going on with the cell parameters. If you run the cell relaxation with Quantum Espresso, you can use my python script “pw2cellvec” to parse all the information you need.

So, first of all, run a “vc-relax” calcuculation with a large (about 2 times bigger) cutoff for both wavefunction (ecutwfc) and charge (ecutrho) and wait for a couple of iterations. These calculations are very expensive, due to the large cutoffs.

Then, just run the script and you will get something like this:

user> pw2cellvec.py -o pwout
- -
AA a= 6.02104498011 b= 6.02109150665 c= 9.45688483193 Vol.= 296.919819153
au a= 11.3781260524 b= 11.3782139748 c= 17.8709224124 c/a= 1.57063846279
alpha= 89.9989981833 beta= 90.0012313781 gamma= 119.99694776
- - Continue reading →

pw2gap: Quick Homo/Lumo gap check for QE

In order to check the Homo/Lumo gap from a Quantum Espresso (pw.x) calculation it is possible to plot the Density Of States and estimate the distance of empty states around the Fermi level (E=0), but what about if we only want a numerical (accurate) value?
This can be done very quick with “pw2gap” script, which can be downloaded from here.

It goes trough all the eigenvalues at the end of the output, including the case of spin-polarization (spin up/down) and k-points, and it simply prints the energies of the HOMO and LUMO orbitals as well as, of course, the difference between them. As simple as that.

It looks like this (for the example of the wustite):

user> pw2gap.py -o pwout
HOMO= 9.4503 LUMO= 11.5898 eV
H/L GAP = 2.1395

sum_states.py: Nice DOS plots from QE outputs

Once, I was making some DOS figures with quantum espresso by using the sumpdos.x program by Andrea Ferretti (included in the QE package), but I felt that I am missing some features. Thereby, I wrote my own script, which does something similar, but it is much more automatized, it is faster to use, and it even produces ready-to-publish graphs automatically.

It is part of the quantum espresso distribution, so if you have the source code of quantum espresso, you have it in the …/*espresso*/PP/tools/sum_states.py. If you don’t have it, or you just want to read the source-code online, you can check for sum_states.py in the quantum espresso repository. Continue reading →

dos-ipr.f: Calculate DOS and IPR with CPMD

The Denstiy Of States (DOS) and the Inverse Partitipation Ratio (IPR) are two interesting properties for understanding the electronic structure of a system.

The DOS is just a histogram counting the amount of states (molecular orbitals/wavefunctions) per energy unit, and analyzing these distributions, we can understand better the electronic behavior of our system.

The IPR is a way to analyze the “amount of localization” of these states, so that the larger the value of IPR the higher the localization around an specific covalent bond.

We used these two properties in an article published in Phys. Rev. B, named “Polymorphism in phase-change materials: melt-quenched and as-deposited amorphous structures in Ge2Sb2Te5 from“, which you can check for more information.  Continue reading →

shrink_traj: Make trajectory files smaller

When we perform molecular dynamics (MD) simulations, sometimes we want to store a frame every time step, in order to increase the quality of the statistics (i.e. when calculating radial distribution functions). But on the other hand, due to their big size, these trajectories might be very difficult to handle by visualization programs (VMD, jmol, …), because they usually load the hole file into memory, and storing every frame (or each few frames) in a long MD with a lot of atoms can produce a trajectory file of several MB or even GB.

The solution: shrinking the big trajectory into a smaller one with my “shrink_traj” script.
This script takes a trajectory file (in xyz format) and copies each nth frame to another file, resulting on a smaller and easier to handle file specially suitable for visualization.