19 September 2013

516. Potential energy diagrammes in latex using pgfplots

Given these recent posts, a post on how to represent reaction coordinates in latex was inevitable.

Note: a much prettier example is found here: https://quantumchemistryniser.wordpress.com/2013/05/25/potential-energy-profile-using-pgfplots/

However, I initially couldn't reproduce it, and after spending a few hours with pgfplots I wanted my own style anyway.


But first, a recap:
Once you've completely left office-like wordprocessors behinds and fully embraced the wonder of typesetting in LaTeX you find yourself wanting to produce more and more of the figures directly in latex, rather than simply importing figures.

While nothing beats gnuplot for plotting, and pov-ray (via e.g. gdis) for fancy figures and bkchem for simple 2d structures, there are some nice and simple latex packages that make life easier: I touch on a few in this post: http://verahill.blogspot.com.au/2012/07/debian-texlive-and-style-files-making.html

For chemical formulae, use mhchem -- I think it's standard in debian for texlive. It's indispensable. Also see http://www.ctan.org/pkg/mhchem

For very, very simple chemical structures you can look at chemfig -- I don't use it much though. http://www.ctan.org/tex-archive/macros/latex/contrib/chemfig/

For drawing simple MO diagrammes, modiagram is very pretty: http://www.ctan.org/pkg/modiagram

For drawing NMR pulse sequences there are a few options: http://verahill.blogspot.com.au/2013/08/499-briefly-drawing-nmr-sequences-using.html and http://verahill.blogspot.com.au/2013/08/498-briefly-drawing-nmr-pulse-sequences.html


A dead end: endiagram
It looks very promising, but it didn't work out for me: http://ctan.mirrorcatalogs.com/macros/latex/contrib/endiagram/

Installation:
mkdir ~/texmf/tex/latex/endiagram -p
cd ~/texmf/tex/latex/endiagram
wget http://ctan.mirrorcatalogs.com/macros/latex/contrib/endiagram/endiagram.sty
sudo texhash

You'll need siunitx and tikz as well, but it seems to come with a standard texlive install on debian. Out of the box I kept getting:

! Undefined control sequence.
<argument> \clist_count:N
\l__endiagram_points_clist

l.16 \ENcurve{1,4,0}

In the end I couldn't get it to work properly, even after installing l3kernel -- I don't think it would've been impossible to figure it out, but I was running out of patience.

The solution: pgfplots
So I kept looking, and found this: https://quantumchemistryniser.wordpress.com/2013/05/25/potential-energy-profile-using-pgfplots/

While I couldn't get it to work 'out of the box' (after deleting the nodes depending on the pdf files), I ended up reading the pfgplots manual: http://pgfplots.sourceforge.net/pgfplots.pdf

pgfplots is not meant for drawing potential energy diagrammes, but it can be made to do it -- it's 'dum' in the sense that it's not aware of what a PE diagram is, so you'll have to draw all objects expliticly. In the end, it's actually quite simple though.

I eventually ended up with these two examples:
The code is here:
\documentclass[10pt,a4paper]{article} \usepackage{pgfplots} \usepackage{tikz} \begin{document} %\begin{center} %\begin{tikzpicture}[scale=1.0] %\draw (0,0) --(1,2); %\end{tikzpicture} %\end{center} \begin{tikzpicture} \begin{axis}[ axis x line= bottom, axis y line= left, xmin=-1, xmax=10, ymin=-1, ymax=6, xlabel=Reaction coordinate, ylabel=$\Delta$E (kcal/mol), xtick=\empty ] \addplot[smooth,solid,blue] coordinates { (0,0)(1,0) (2,4)(3,4) (4,3)(5,3) (6,4)(7,4) (8,1)(9,1)}; \node at (axis cs:4.5,2.5) {Intermediate}; \end{axis} \end{tikzpicture} \begin{tikzpicture} \begin{axis}[ axis x line= bottom, axis y line= left, xmin=-1, xmax=6, ymin=-1, ymax=15, xlabel=Reaction coordinate, ylabel=$\Delta$E (kcal/mol), xtick=\empty ] \addplot[smooth,solid,blue] coordinates { (0,0)(1,0)}; \addplot[smooth,dotted,black] coordinates {(1,0)(2,11.7)}; \addplot[smooth,solid,blue] coordinates {(2,11.7)(3,11.7)}; \addplot[smooth,dotted,black] coordinates {(3,11.7)(4,10.5)}; \addplot[smooth,solid,blue] coordinates {(4,10.5)(5,10.5)}; \end{axis} \end{tikzpicture} \end{document}


17 September 2013

515. Very briefly: boot issues and nvram (anecdote)

I had to move offices at work last Friday (which partly explains my radio silence).

One of my nodes, with a Biostar 880G+ motherboard, wouldn't boot. The issue manifested itself in a few different ways (and not always at the same time):

Most of the time the bios portion of the boot would be very slow, and I'd get a loooooong beep when the voltages were echoed during what I presume is POST. I never made it to the GRUB screen. I'd had this happen once or twice before, especially when the computer had been turned off and/or unplugged for a longer period of time.

Occasionally I'd get a message about NVRAM being checked..

After
1) having checked that nothing had vibrated loose during the transport
2) having removed all cards and drives except a single RAM stick
3) having checked that RAM stick
I still wasn't getting any closer to getting booted. So I went into the BIOS. Having seen the message about nvram and having done a bit of reading online I selected to erase the NVRAM.

And suddenly booting worked...once.

After that I kept getting messages about nvram being corrupted and boot fail, or 'NVRAM...Update failed'.

I then set the CMOS jumpers to pin 2 and 3 (small jumper next to the CMOS battery)

And then it worked! It booted just fine! But I kept getting the following error echoing on my screen:

[ 120.718683] hub 3-0:1.0: unable to enumerate USB device on port 2 [ 120.932378] usb 3-2: new full-speed USB device number 28 using ohci_hcd [ 121.072497] usb 3-2: device descriptor read/64, error -62 [ 121.316673] usb 3-2: device descriptor read/64, error -62 [ 121.556904] usb 3-2: new full-speed USB device number 29 using ohci_hcd [ 121.697019] usb 3-2: device descriptor read/64, error -62 [ 121.941224] usb 3-2: device descriptor read/64, error -62 [ 122.181425] usb 3-2: new full-speed USB device number 30 using ohci_hcd
Unplugged all the USB devices. Turned out to be the KVM (ritmo). Turning it off (i.e. unplugging all USB input to the KVM so that it lost power) and plugging it in again solved it.

And ever since then all has been fine.

Presumably it's time to change the CMOS battery.

514. Extracting Frequency data from a gaussian 09 calculation for gnuplot

This is another python script.

Say you've done a computation along the lines of this:
#P rBP86/GEN 5D Pseudo(Read) Opt=() Freq=() SCF=(MaxCycle=999 ) Punch=(MO) Pop=()
and want the data in a neat data file, like this:
33.237 0.0023 0.0536 39.9976 0.0043 0.8305 69.7345 0.0129 0.3348 84.7005 0.0173 0.7027 [..] 3133.0068 6.2938 0.6114 3143.8021 6.3551 0.3775 3164.9242 6.4685 0.8829 3221.8787 6.6972 4.6005

Then you can use the following python (2.x) script, g09freq:

#!/usr/bin/python
# Compatible with python 2.7 
# Reads frequency output from a g09 (gaussian) calculation
# Usage ex.: g09freq g09.log ir.dat
import sys 

def ints2float(integerlist):
    for n in range(0,len(integerlist)):
        integerlist[n]=float(integerlist[n])
    return integerlist

def parse_in(infile):
    g09output=open(infile,'r')
    captured=[]
    for line in g09output:
        if ('Frequencies' in line) or ('Frc consts' in line) or ('IR Inten' in line):
            captured+=[line.strip('\n')]
    g09output.close()
    return captured
    
def format_captured(captured):
    vibmatrix=[]
    steps=len(captured)
    for n in range(0,steps,3):
        freqs=ints2float(filter(None,captured[n].split(' '))[2:5])
        forces=ints2float(filter(None,captured[n+1].split(' '))[3:6])
        intensities=ints2float(filter(None,captured[n+2].split(' '))[3:6])
        for m in range(0,3):
            vibmatrix+=[[freqs[m],forces[m],intensities[m]]]
    return vibmatrix

def write_matrix(vibmatrix,outfile):
    f=open(outfile,'w')
    for n in range(0,len(vibmatrix)):
        item=vibmatrix[n]
        f.write(str(item[0])+'\t'+str(item[1])+'\t'+str(item[2])+'\n')
    f.close()
    return 0

if __name__ == "__main__":
    infile=sys.argv[1]
    outfile=sys.argv[2]

    captured=parse_in(infile)

    if len(captured)%3==0:
        vibmatrix=format_captured(captured)
    else:
        print 'Number of elements not divisible by 3 (freq+force+intens=3)'
        exit()
    success=write_matrix(vibmatrix,outfile)
    if success==0:
        print 'Read %s, parsed it, and wrote %s'%(infile,outfile)


Run it as e.g.
g09freq g09.log test.out

The output is compatible with gnuplot:
gnuplot> set xrange [3500:0]
gnuplot> set yrange [10:-1]
gnuplot> plot './test.out' u 1:2 w impulse



It's trivial to add gaussian broadening (see e.g. this post)