OS-X Mavericks for Scientific Computing

My Macbook Pro (Early 2011) running Snow Leopard was starting to get a little slow, so I decided to upgrade it with a 250Gb SSD, and use the opportunity of an empty disk to give a chance to OS-X Mavericks.
In this post, I am going to explain step by step how did I set up my laptop with Mavericks and some extra tools for enabling scientific computing and computational chemistry tasks.
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

agr2gpt: convert xmgrace to gnuplot or mathplotlib

Some people seem to think that xmgrace is somehow cool, maybe because it contains an “x” on its name. Y personally prefer gnuplot or mathplotlib, but I got recently some “.agr” files I had to edit, and I finally decided to write a python script called agr2gpt, which makes a simple conversion from xmgrace files to gnuplot script, and extracts the data to different files (one per curve).

This is just the first version, so it just handles basic 2D graphs with several curves.
Besides the x and y data for each curve (unlimited amount of curves), it also parses the “xlabel” the “ylabel” and the legend of each curve.

As in most of my scripts, the “-h” option will provide some extra information about the usage, but in this first version, there are only two options: “-f” for the file name, and “-pl” if you want the pylab window to pop up with the image. 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

Automate molecular figures

One of the things I like most from being a computational chemist, is the possibility to automate tasks by using scripts.
There is usually no problem for automating almost anything, but unfortunately, there are a few tricky to automate tasks, for example, when it comes to generate atomic/molecular figures.
Let’s think that we have to calculate the energies of a lot of different molecules, and after doing that, we certainly need to check the results. We can easily “grep” the total energy (or any other numerical/ASCII value) from the output files but usually, the only way to take a look to the geometry of the structures is by going trough all the outputs and opening them with a graphic molecular viewing program (such as molden, vmd, xmakemol, jmol, …). This task can be quite painful, specially for the geek people who does not like wasting their energy using the mouse if it is not completely necessary. Continue reading

now for different supercomputers

I already shared my now script intended for a small “in home” computing cluster. This time, I am sharing two other variations of that script, designed for the queuing systems of two different supercomputers: Mare Nostrum at BSC, Juropa at JSC and HLRN.

Both of them are the predecessors of the newer now script, and they only have the very basic features, but if you need an script that works right out of the box for any of the supercomputers mentioned above, this is the quick solution.

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 OStates (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