Blog: Experiment: Conformational analysis of ethane: C2H6_rotate-bond.py

File C2H6_rotate-bond.py, 2.1 KB (added by FrederikHeber, 8 years ago)

python script using pyMoleCuilder for rotating ethane around its C-C bond

Line 
1import pyMoleCuilder
2import os, sys, csv
3
4# check arguments
5if len(sys.argv) < 6:
6 sys.stderr.write("Usage: "+sys.argv[0]+" <source file> <steps> <server> <port> <workdir> <energyfile> <conformationfile>\n")
7 sys.exit(255)
8
9sourcefile=sys.argv[1]
10steps=int(sys.argv[2])
11serveraddress=sys.argv[3]
12serverport=sys.argv[4]
13workdir=sys.argv[5]
14energyfile=sys.argv[6]
15conformationfile=""
16if len(sys.argv) >= 7:
17 conformationfile=sys.argv[7]
18
19# ========================== Stored Session BEGIN ==========================
20pyMoleCuilder.MoleculeLoad(sourcefile)
21pyMoleCuilder.ParserSetParserParameters("mpqc", "basis=6-31G;theory=MBPT2;")
22energies={}
23stepwidth=int(360/steps)
24for step in range(0,360,stepwidth):
25 pyMoleCuilder.wait()
26 # select bond
27 pyMoleCuilder.SelectionClearAllAtoms()
28 pyMoleCuilder.SelectionAtomByElement("C")
29 # advance to next time step
30 pyMoleCuilder.WorldStepWorldTime("1")
31 # rotate bond
32 pyMoleCuilder.MoleculeRotateAroundBond(str(stepwidth), "1")
33 # calculate energy
34 pyMoleCuilder.SelectionAllAtoms()
35 pyMoleCuilder.FragmentationFragmentation("BondFragment", "3", "2", "1", "1", "", "0", "5", "0", "0", "0")
36 pyMoleCuilder.FragmentationFragmentationAutomation(serveraddress, serverport, "", "0", "5", "3", "3", "0", "0", "0", "0", "0")
37 pyMoleCuilder.FragmentationAnalyseFragmentationResults("0", "BondFragment", "0")
38 # extract energy and place in hash map
39 pyMoleCuilder.wait()
40 with open(workdir+"/BondFragment_Energy.dat", 'rb') as csvfile:
41 energyreader = csv.reader(csvfile, delimiter='\t', quotechar='|')
42 for row in energyreader:
43 if (row[0]=="2"):
44 energies[step]=row[1]
45 break
46# output all
47if len(conformationfile) != 0:
48 pyMoleCuilder.WorldOutputAs(conformationfile)
49# save energies as csv
50with open(energyfile, 'wb') as csvfile:
51 energywriter = csv.writer(csvfile, delimiter=',', quotechar='|', quoting=csv.QUOTE_MINIMAL)
52 energywriter.writerow(["step", "energy"])
53 for key in sorted(energies.keys()):
54 energywriter.writerow([str(key), str(energies[key])])
55# =========================== Stored Session END ===========================
56
57sys.exit(0)