Adventures in Computational Chemistry

Recently I have been working on Linear combination of Distance (LCOD). One has to extract a set of structures (may vary from 100 to 1000 to any bigger number depending on requirement, accuracy) from a very long trajectory. If one is looking for very specific orientation, then it is necessary to go through the whole trajectory and choose the required frames and select that geometry, which is very tedious and time taking job. But more often a random selection from the trajectory is done to take care of all the diversity in the trajectory. In Amber, one can use ptraj to select a frame from trajectory and created its restart file to launch LCOD calculation. To launch such 100 or 1000 calculation one need to modify or create 100 or 1000 such input files to feed to ptraj with the frame number and the name of the restart file. I am using a dictionary in python to do so. Beginning with 2 lists: one having frame numbers to be selected another with name of the restart files correspondin

I am still struggling with the after effects of Linux reinstallation on my system. Two major problems: 1) Firefox and openoffice is still not in English. Its giving me hard time. 2) The functions in openoffice  version of excel is not working. I am yet to sort that out. In the mean time, I desperately needed it to do simple task. subtracting values of one column from one in another column and write the difference in next column. So, I had input in this form: ------------------ 0    2.179139    0.951816 1    2.314531    1.006480 2    2.191692    0.951938 3    2.079704    0.952747 4    2.092761    0.939545 5    2.362056    0.984791 6    2.153595    0.957934 7    2.213173    0.951768 8    2.175687    0.953091 -------------------- I needed to subtract  column 3 data from column 2 data. So, I tried this: --------------------------- #!/usr/bin/env python fo=open('name.txt' +'_output','w') with open('name.txt') as fd:  for line i

For some reason our administrator decided to update OS of my desktop. He installed latest version Ubuntu 11.10. My first impression was not very good. I did not find it userfriendly. Main problem is accessing different windows you open. Say, you open a firefox window with many work related tabs opned in it. This window will minimize to the left sidebar under firefox icon that also is the desktop icon to access firefox. Now, if you want open another firefox window which should have pages not work related but say personal like gmail or youtube, you click on the same icon. Instead of what you wish, the same firefox window with work related pages which you minimized few moments back reappears. So, to distinguish work related pages from personal pages is not so stright forward. One option is you open your personal page in one of the tabs in work related firefox window and then drag and drop that tab out of the window. Then it allows you to have another firefox window. But why is needs

GNUPLOT can be used to generate histograms. The script below allows one to plot variations in angle from 0° to 180°. Each histogram box width correspond to 3° (covered here under bw option). To execute this script; type this in the terminal: gnuplot histogram.p You should be in the folder where below script is saved in a file name "histogram.p" and a dat file (here, it is N---H-O.dat) should be present. histogram.p ------------------ set term gif set output "angle_N---H-O.gif" set xrange [0.0:180.0] set boxwidth 2.0 absolute set style fill solid 1.0 border -1 set xtics 10 set title "molecule 3.0A Angle N---H-O" set ylabel "occurance" set xlabel "angle" bw = 3   # substitute what you want bin(x,width)=width*floor(x/width) plot 'N---H-O.dat' using (bin($2,bw)):(1.0) smooth freq with boxes ------------------ The dat file is an input file with data in the format: 0    23.680344 1    73.766251 2    67.4

I needed to delete all the lines in a set of files except for for the first forty. Since the number of files was big, I used this script in python to do the job: -------------------------------- !/usr/bin/python #script to read files and delete lines import os listing=os.listdir(" /home/kjoshi/Documents/PY/useful/file_edit ") for infile in listing:  print infile  if infile == "files_line_delete.py":    print "This is NOT a good file"  else:   fi=open(infile).readlines()   fo=open(infile +'_output','w')   del fi[40: ]   fo.writelines(fi) --------------------------------- Task is simple. Keep all the files needed to be edited in one folder along with this script. Substitute the path shown in bold with proper path of your folder. Using command $python files_line_delete.py you get output files with each file having first forty lines.

The task was to create a set of folders say 1,2,3,4 and in each of these four folders create four sub folders. Thus 16 folders were to be created and in each of these 16 folders a common set of files were to be copied. I used c shell scripting to do the job. I needed two lists. One list for primary folders list1: 1 2 3 4 and another list for the subfolders. list2: 5 6 7 8 It worked as: ------------------- bash-3.2$ csh r410comp2% r410comp2%foreach f (`cat ../list1`) foreach? echo $f foreach? mkdir "$f" foreach? cd "$f" foreach? foreach g (`cat ../../list2`) foreach? echo $g foreach? mkdir "$f"_"$g" foreach? cd "$f"_"$g" foreach? cp /home/kjoshi/W/work/complex/2_ion/umbrella/6/6_75/u_6_75.in . foreach? cp /home/kjoshi/W/work/complex/2_ion/umbrella/6/6_75/amber11_intel_ref_md . foreach?cp /home/kjoshi/W/work/complex/2_ion/umbrella/6/6_75/2_ion.inpcrd . foreach?cp /home/kjoshi/W/work/complex/2_io

This script allows one to exchange two columns in a file.And it do so for all the files in a directory. The python script is : ---------------------------- #script to read and execute column exchange for all files in the directory import os listing=os.listdir(" /home/K/Documents/PY/column_exchange ") for infile in listing:  print infile  if infile !="column_interchange_script.py":   fi=open(infile).readlines()   fo=open(infile +'_output','w')   space="   "   for line in fi:    first,second=line.split()    fo.write('%s%s%s\n'%(second,space,first))  else:    print "This is not interesting file"     ----------------------------------------------------    Some expert might simplify it a lot more, but for new python user like me, this is a big leap. So, how does it work: Create a directory where in all the files for which column exchange need to be done are stored. Keep this python script also in the sa

Now that we have the lib file with the information for all the molecules we need to include and also the corresponding combined pdb file, creating prmtop and inpcrd files should not be a difficult task. So, to begin with we have a lib file name: molecule.lib a pdb file for 4 molecules name: four.pdb The  task is to generate a prmtop and inpcrd file. Initiate xleap using: $AMBERHOME/exe/xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99 In xleap window > source leaprc.gaff # load the lib file by: > loadoff molecule.lib I am not sure but some times it helps to also load the frcmod files we created for the four molecules the information of which is stored in the lib file. > loadamberparams aldehyde.frcmod # Load the pdb file: > P=loadpdb four.pdb You should see some thing like this in the xleap window now: ------------------- > P=loadpdb four.pdb Loading PDB file: ./four.pdb Warning: name change in pdb file residue 1 ;  this residue is split i

Packmol is a software that allows you to generate a single combined PDB file for all the molecules you are considering. It is a free software and can be downloaded from this link: http://www.ime.unicamp.br/~martinez/packmol/ The installation guide is also provided on the page and its really easy to follow. http://www.ime.unicamp.br/~martinez/packmol/ You can also find some examples but I tried on my own. Anyway what you need is an input file which looks like: ------------------ tolerance 2.0 output trial.pdb add_amber_ter filetype pdb structure 2a.pdb   number 1   inside cube 0. 0. 0. 10. end structure structure 2d.pdb   number 1   inside cube 0. 0. 0. 10. end structure --------------------------------- IT includes some information like: how close the molecules considered can be.(tolerance) name of output pdb file "add_amber_ter" adds term "TER" after every molecule thus defining the start and end of every molecule. After mentioning the

Sometimes you totally fail to go ahead with the prepin and frcmod files you have generated from the newly generated PDB file even after you follow the procedure given in previous two posts. And you keep wondering? In such a case you need a more detailed procedure. First take a look at the pdb file. For example: Here is one such pdb file generated from the gabedit that troubled me for some time. --------------- HEADER    PROTEIN COMPND    UNNAMED AUTHOR    GENERATED BY GABEDIT 2.3.8 at Tue Jul 26 16:53:48 2011 ATOM  1      C   C     1       -2.907  -0.997  -0.385  1.00  300.00          C 0.0000   ATOM  2      N   N     1       -2.060  -0.981  -1.189  1.00  300.00          N 0.0000   ATOM  3      C   C     1       -1.026  -0.973  -2.172  1.00  300.00          C 0.0000   ATOM  4      H   H     1       -0.409  -1.865  -2.044  1.00  300.00          H 0.0000   ATOM  5      H   H     1       -0.409  -0.082  -2.038  1.00  300.00          H 0.0000   ATOM  6      H   H     1      

In earlier post we saw how we create pdb files for the structures we want to consider for multi component molecular dynamics analysis. So now we have say PDB files for four structures - an aldehyde, nitrile, amine and a phenol molecule. Now, In the first step we generate the prepin files for all these molecules by  using: $AMBERHOME/exe/antechamber -i sustiva.pdb -fi pdb -o sustiva.prepin -fo prepi -c bcc -s 2 Here the molecules are neutral in nature. If you are working with ionic systems you need to add -nc flag in the end followed by charge of the moelcule. So, if the moelcule has +1 charge, the command will be: $AMBERHOME/exe/antechamber -i sustiva.pdb -fi pdb -o sustiva.prepin -fo prepi -c bcc -s 2 -nc 1 If it is -1, command is: $AMBERHOME/exe/antechamber -i sustiva.pdb -fi pdb -o sustiva.prepin -fo prepi -c bcc -s 2 -nc -1 If you have followed earlier preocedure properly then your prepin file should be created in no time. You can see some lines of a prepin file as

If Gabedit is installed and functioning on your system, type gabedit  in terminal. Gabedit will open and you can find an option "Geometry" in the tool bar. Select "Geometry"--> "Draw" After some trials and playing you will finally learn how to generate a structure. Once your structure is ready, you need to set the atom types and charges. For this click on "M" button in the top left of the "Draw Geometry" window and select options as: M-->"Set"--> "Atom Type & charge using PDB template" Here if you use "Atom Types using connection types", it helps. You can check if the atom type and charge has been set according to the amber requirements by an indirect method. Try to run a quick MM calculation. Run MM calculations by clicking: "M"-->Molecular Mechanics-->optimization It may sometimes end up with a pop up window that suggests that either atom type or charge is not

I have a system with four molecules and the initial challange is to generate a combined PDB file for the four molecules together which is compatible with the AMBER software. Apparently it is a simple problem.  But here it is a problem because pdb version of molecules generated using different build up softwares find compatibility issue with Amber. Also I found very rare examples where people are trying to study multicomponent system using Amber. In other words even the Amber mailing list (http://archive.ambermd.org/) was not very helpful in solving the problem. I tried many different possibilities. I used different tutorials for Amber, but none of them were direct solutions. I finally created my own procedure understanding these tutorials and combining them together. In coming few I will try to explain how you can generate amber compatible PDB file for multi component system

In this new addition I will be working on totally new problems using totally new approaches and totally new softwares. My work will be revolving around Organic reaction mechanism I will be doing molecular dynamics I should mainly be using AMBER and I don't know for some reason may be Gaussian, CPMD etc. I am going to use a new interface called "gambit" to generate geometries and the input files. Initial task is to learn Amber. Got ambertools installed and now I have to go through the tutorials.