11 October 2013

519. Formatting an XYZ molecular geometry file using python

This is a silly little script -- ECCE, which I use to manage all my computations, is very particular about what XYZ files it can and cannot read -- if the number of spaces between the coordinates are wrong, it will fail to read the structure. So here's a python script, based on parts I've cannibalised from other scripts that I've written in the past (i.e. it could be made more elegant, but I had the parts ready and just wanted something that works) which takes a somewhat ugly XYZ file and turns it into something beautiful. At least in the eyes of ECCE.

The impetus for this was that someone gave me an XYZ file the had generated by copying columns from Excel. On Mac. mac2unix took care of the line endings, but there were tabs (\t) all over the place.

Usage:
polish_xyz ugly.xyz pretty.xyz

Script:
#!/usr/bin/python
import sys 

def getrawdata(infile):
    f=open(infile,'r')
    n=0 
    preamble=[]
    struct=[]
    
    for line in f:
        if n<2: data-blogger-escaped-if="" data-blogger-escaped-line.rstrip="" data-blogger-escaped-n="" data-blogger-escaped-preamble="">3:
            line=line.rstrip()
            struct+=[line]
        n+=1
    xyz=[struct]
    return xyz, preamble

def genxyzstring(coords,elementnumber):
    x_str='%10.5f'% coords[0]
    y_str='%10.5f'% coords[1]
    z_str='%10.5f'% coords[2]
    element=elementnumber
    xyz_string=element+(3-len(element))*' '+10*' '+\ 
    (8-len(x_str))*' '+x_str+10*' '+(8-len(y_str))*' '+y_str+10*' '+(8-len(z_str))*' '+z_str+'\n'
    
    return xyz_string

def getstructures(rawdata,preamble,outfile):
    n=0 
    for structure in rawdata:
        n=n+1
        num="%03d" % (n,)
        g=open(outfile,'w')
        itson=False
        cartesian=[]
    
        for item in structure:
            coordx=filter(None,item.split(' ')) 
            coordy=filter(None,item.split('\t'))
            if len(coordx)>len(coordy):
                coords=coordx
            else:
                coords=coordy
            coordinates=[float(coords[1]),float(coords[2]),float(coords[3])]
            element=coords[0]
            cartesian+=[genxyzstring(coordinates,element)]
        g.write(str(preamble[0])+'\n')
        g.write(str(preamble[1])+'\n')
        for line in cartesian:
            g.write(line)
        g.close()
        cartesian=[]
    return 0

if __name__ == "__main__":
    infile=sys.argv[1]
    outfile=sys.argv[2]
    xyz,preamble=getrawdata(infile)
    structures=getstructures(xyz,preamble,outfile)

518. Generating enantiomers of molecular structures given in XYZ coordinates using python.

What I'm showing here is probably overkill -- there may be better ways of doing this with sed/awk*. However, since I had most of the code below ready as it was part of another script, doing this is python was quick and easy. Plus it's portable-ish.

*[ECCE, which I use for managing my calculations, is very, very particular about the format of the XYZ file, including the number of chars between coordinates. So simply doing e.g.

cat molecule.xyz|gawk '{print $1,$2,$3,-$4}'

won't work. Not all pieces of software are that picky when it comes to xyz coordinates though -- jmol, for example, is very forgiving.]

Anyway, the script below flips a molecular structure by taking the geometry given as a (properly formatted) XYZ file, and changing the sign in front of the Z coordinates. It's that simple.
Save it, call it flip_xyz, make it executable and call it using e.g.

flip_xyz molecule.xyz flipped_molecule.xyz

Script:
#!/usr/bin/python
import sys 

def getrawdata(infile):
    f=open(infile,'r')
    n=0 
    preamble=[]
    struct=[]
  
    for line in f:
        if n<2: data-blogger-escaped-if="" data-blogger-escaped-line.rstrip="" data-blogger-escaped-n="" data-blogger-escaped-preamble="">1:
            line=line.rstrip()
            struct+=[line]
        n+=1
    xyz=[struct]
    
    return xyz, preamble

def genxyzstring(coords,elementnumber):
    x_str='%10.5f'% coords[0]
    y_str='%10.5f'% coords[1]
    z_str='%10.5f'% -coords[2]
    element=elementnumber
    xyz_string=element+(3-len(element))*' '+10*' '+\ 
    (8-len(x_str))*' '+x_str+10*' '+(8-len(y_str))*' '+y_str+10*' '+(8-len(z_str))*' '+z_str+'\n'
    
    return xyz_string

def getstructures(rawdata,preamble,outfile):
    n=0 
    for structure in rawdata:
        n=n+1
        num="%03d" % (n,)
        g=open(outfile,'w')
        cartesian=[]
    
        for item in structure:
            coordx=filter(None,item.split(' ')) 
            coordy=filter(None,item.split('\t'))
            if len(coordx)>len(coordy):
                coords=coordx
            else:
                coords=coordy
    
            coordinates=[float(coords[1]),float(coords[2]),float(coords[3])]
            element=coords[0]
            cartesian+=[genxyzstring(coordinates,element)]
    
        g.write(str(preamble[0])+'\n')
        g.write(str(preamble[1])+'\n')
        for line in cartesian:
            g.write(line)
        g.close()
        cartesian=[]
    return 0

if __name__ == "__main__":
    infile=sys.argv[1]
    outfile=sys.argv[2]
    xyz,preamble=getrawdata(infile)
    structures=getstructures(xyz,preamble,outfile)


20 September 2013

517. Very briefly: Prime95 (GIMPS) on linux

I'm very unhappy about a newly built node which randomly crashes and reboots when running long jobs. More about that later, but here are the specs: FX 8350, 4x8 Gb RAM GSkill Ripjaws, ASRock FX990 Extreme3, Corsair GS700, MSI N210, ASUS NX1101 in an Antec GX700 case, running Wheezy with stock kernel (3.2.0-4 amd64).

I've tested the RAM using memtest86+  and found no errors, the rig uses a 700 W Corsair PSU which /should/ provide enough power, and I see no evidence of overheating based on a cronjob which runs every 2 minutes. Anyway, the first step in troubleshooting is finding a good way of reproducing the error reliably, and prime95 is what the windows overclockers use to stresstest.

Turns out prime95 (actually GIMPS) can run in a few different modes which tests different aspects of you system, which makes it sound like a pretty good program for my purposes.

See here for more information: http://www.mersenne.org/freesoft/

mkdir ~/tmp/mprime -p
cd ~/tmp/mprime
wget http://www.mersenne.info/gimps/p95v279.linux64.tar.gz
tar xvf p95v279.linux64.tar.gz
./mprime
Welcome to GIMPS, the hunt for huge prime numbers. You will be asked a few simple questions and then the program will contact the primenet server to get some work for your computer. Good luck! Attention OVERCLOCKERS!! Mprime has gained a reputation as a useful stress testing tool for people that enjoy pushing their hardware to the limit. You are more than welcome to use this software for that purpose. Please select the stress testing choice below to avoid interfering with the PrimeNet server. Use the Options/Torture Test menu choice for your stress tests. Also, read the stress.txt file. If you want to both join GIMPS and run stress tests, then Join GIMPS and answer the questions. After the server gets some work for you, stop mprime, then run mprime -m and choose Options/Torture Test. Join Gimps? (Y=Yes, N=Just stress testing) (Y): N Number of torture test threads to run (3): 2 Choose a type of torture test to run. 1 = Small FFTs (maximum FPU stress, data fits in L2 cache, RAM not tested much). 2 = In-place large FFTs (maximum heat and power consumption, some RAM tested). 3 = Blend (tests some of everything, lots of RAM tested). 11,12,13 = Allows you to fine tune the above three selections. Blend is the default. NOTE: if you fail the blend test, but can pass the small FFT test then your problem is likely bad memory or a bad memory controller. Type of torture test to run (3): 1 Accept the answers above? (Y): Y [Main thread Sep 20 11:06] Starting workers. [Worker #1 Sep 20 11:06] Worker starting [Worker #1 Sep 20 11:06] Setting affinity to run worker on any logical CPU. [Worker #2 Sep 20 11:06] Worker starting [Worker #2 Sep 20 11:06] Setting affinity to run worker on any logical CPU. [Worker #1 Sep 20 11:06] Beginning a continuous self-test to check your computer. [Worker #1 Sep 20 11:06] Please read stress.txt. Hit ^C to end this test. [Worker #2 Sep 20 11:06] Beginning a continuous self-test to check your computer. [Worker #2 Sep 20 11:06] Please read stress.txt. Hit ^C to end this test. [Worker #1 Sep 20 11:06] Test 1, 180000 Lucas-Lehmer iterations of M580673 using AMD K10 type-1 FFT length 28K, Pass1=112, Pass2=256. [Worker #2 Sep 20 11:06] Test 1, 180000 Lucas-Lehmer iterations of M580673 using AMD K10 type-1 FFT length 28K, Pass1=112, Pass2=256. CTRL+C
And so on.