Jump to content

Talk:Mod:Hunt Research Group/build freq file

From ChemWiki

build a frequency input file

  • this script will need editing for your particular job
  • the script extracts the last geometry from the log file supplied and builds a com file, useful especially if you are using gen and pseudo-potentials
  • copy the contents below into a file "build_freq.py"
  • run with "python build_freq.py file_name" give the logfile without the .log extension
#
# build_freq.py 
#
# python script to extract the final structure from a gaussian log file and build a frequency input file
# to run the script type python build_freq.py input_file_name (without extension)
#
import sys
import os
dir=os.getcwd()
#
if len(sys.argv) == 1:
  print('to run the script please type: python build_freq.py input_file_name \(without extension\)')
  sys.exit()
else:
  base=str(sys.argv[1])
  file=base+'.log'
  out=base+'.xyz'
  base2=base.replace('opt','freq')
  com=base2+'.com'
  s='directory is '
  print('{0:}{1:}'.format(s,dir))
  s='input file is '
  print('{0:}{1:}'.format(s,file))
  s='output file is '
  print('{0:}{1:}'.format(s,com))
#close if
#
# set some defaults
col_num=[]
col_atom=[]
col_type=[]
col_x=[]
col_y=[]
col_z=[]
rep=0
go=False

# open file and read lines one by one
f=open(file,"r")

line =f.readline()
while line :
  line =f.readline()

#==extract GEOMETRY ==
# check for standard orientation text set go if the text appears
# this is the first geometry 
# rep is number of steps ie geometries
# count will be the number of atoms
# n is stepping through the lines in std orientation section
  if 'Standard orientation:' in line:
    go=True
    s='Structure in standard orientation found '
    if rep == 0 : print('{0:}'.format(s))
    n=1
    rep=rep+1
    count=0
# close if
#
# if go is on then store the structure
# skip the next 4 lines, on the 5th start collecting structure
# don't print the lines of dashes
# only print if you have not got to rot const line
  while go:
    n = n+1
    m = 0
    line =f.readline()
    if 'Rotational constants (GHZ):' in line: m=-1    
    if '---------' in line: m=1
    if n > 4 and m == 0 :
#      print('line: {0:}'.format(line))
      count=count+1
      a1,b1,c1,d1,e1,f1=line.rstrip().split()
      col_num.append(int(a1))
      col_atom.append(int(b1))
      col_type.append(int(c1))
      col_x.append(float(d1))
      col_y.append(float(e1))
      col_z.append(float(f1))
    if m == -1 : go=False
# close while
#close while

#now build com file
#
#first give all the comp details
#
c=open(com,"w")
eol="\n"
s="%nprocshared=40\n%mem=120000MB\n"
c.write('{0:}'.format(s))
base2=base.replace('opt','freq')
s='%chk='+base2+'.chk'
c.write('{0:}{1:}'.format(s,eol))
s="# freq b3lyp/gen pseudo=cards empiricaldispersion=gd3bj geom=cartesian\n gfinput int=ultrafine scf=conver=9 test\n\n"
c.write('{0:}'.format(s))
s="Title Card Required\n\n0 1\n"
c.write('{0:}'.format(s))
#
# now put in the coordinates
# we may have muliple geometries print only the last one
#
# n is used to count backward from last entry to find start atom of last geometry
# count will be the number of atoms
# n then becomes the position marker in the geometry array
# count then becomes the final position marker in the geometry array
# rep is number of steps ie geometries
n=0
if rep != 1 :
  n=len(col_num)-count
  count=count*rep
  s=' structures found only the last is saved'
  print('{0:}{1:}'.format(rep,s))
# endif
while n < count:
  c.write('{0:>4d}  {1: 09.6f}  {2: 09.6f}  {3: 09.6f} \n'.format(col_atom[n],col_x[n],col_y[n],col_z[n]))
#  print('{0:>4d} {1:>4d}  {2: 09.6f}  {3: 09.6f}  {4: 09.6f}'.format(col_num[n],col_atom[n],col_x[n],col_y[n],col_z[n]))
  n=n + 1
#close while

# now put in the gen information
s="\nC N H Cl 0\n6-311+G(d,p)\n****\n"
t="Sn 0\nLANL2DZ\n****\n\n"
u="Sn 0\nLANL2DZ\n\n\n"
c.write('{0:}{1:}{2:}'.format(s,t,u))

# close all open files
f.close()
c.close()
#