Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Modify molecule.py to support Q-Chem EDA/FDA calculations. #192

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
37 changes: 28 additions & 9 deletions geometric/molecule.py
Original file line number Diff line number Diff line change
Expand Up @@ -1381,8 +1381,8 @@ def __deepcopy__(self, memo):
# (Q-Chem input section name) and the second element
# is a list (Q-Chem input section data).
New.Data[key] = []
for SectName, SectData in self.Data[key]:
New.Data[key].append((SectName, SectData[:]))
for SectName, SectData, SectInsert in self.Data[key]:
New.Data[key].append((SectName, SectData[:], SectInsert))
else:
raise RuntimeError("Failed to copy key %s" % key)
return New
Expand Down Expand Up @@ -3595,11 +3595,14 @@ def read_qcin(self, fnm, **kwargs):
mult = 0
Answer = {}
SectionData = []
SectionInsert = {}
template_cut = 0
section_ln = 0 # Line number within a section
readsuf = True
suffix = [] # The suffix, which comes after every atom line in the $molecule section, is for determining the MM atom type and topology.
ghost = [] # If the element in the $molecule section is preceded by an '@' sign, it's a ghost atom for counterpoise calculations.
infsm = False
found_charge_mult = False

# The 'iso-8859-1' prevents some strange errors that show up when reading the Archival summary line
for line in open(fnm, encoding='iso-8859-1').readlines():
Expand All @@ -3610,6 +3613,7 @@ def read_qcin(self, fnm, **kwargs):
zmatrix = True
if re.match(r'^\$',line):
wrd = re.sub(r'\$','',line).lower()
section_ln = 0
if zmatrix:
if wrd == 'end':
zmatrix = False
Expand All @@ -3630,7 +3634,7 @@ def read_qcin(self, fnm, **kwargs):
qcrem = OrderedDict()
if reading_template:
if section != 'external_charges': # Ignore the external charges section because it varies from frame to frame.
template.append((section,SectionData))
template.append((section,SectionData,SectionInsert))
SectionData = []
else:
section = wrd
Expand Down Expand Up @@ -3658,14 +3662,15 @@ def read_qcin(self, fnm, **kwargs):
if readsuf and len(sline) > 4:
whites = re.split('[^ ]+',line)
suffix.append(''.join([whites[j]+sline[j] for j in range(4,len(sline))]))
elif re.match("[+-]?[0-9]+ +[0-9]+$",line.split('!')[0].strip()):
if not found_first_mol:
charge = int(sline[0])
mult = int(sline[1])
elif re.match("[+-]?[0-9]+ +[0-9]+$",line.split('!')[0].strip()) and not found_charge_mult:
charge = int(sline[0])
mult = int(sline[1])
found_charge_mult = True
elif infsm and re.sub(r'\$','',line).lower() != 'end':
pass
else:
SectionData.append(line)
SectionInsert[section_ln] = line
# SectionData.append(line)
elif reading_template:
if section == 'basis':
SectionData.append(line.split('!')[0])
Expand All @@ -3677,6 +3682,7 @@ def read_qcin(self, fnm, **kwargs):
qcrem[S[0].lower()] = ''.join(S[2:])
else:
SectionData.append(line)
section_ln += 1
elif re.match('^@+$', line) and reading_template:
template.append(('@@@', []))
elif re.match('Welcome to Q-Chem', line) and reading_template and found_first_mol:
Expand Down Expand Up @@ -4209,6 +4215,14 @@ def write_qcin(self, **kwargs):
read = kwargs['read']
else:
read = False

def insert_section(linelist, section_ln, sect_insert):
section_ln += 1
while section_ln in sect_insert:
linelist.append(sect_insert[section_ln])
section_ln += 1
return section_ln

for SI, I in enumerate(selection):
fsm = False
remidx = 0
Expand All @@ -4220,7 +4234,7 @@ def write_qcin(self, **kwargs):
for i in range(len(extchg)):
out.append("% 15.10f % 15.10f % 15.10f %15.10f" % (extchg[i,0],extchg[i,1],extchg[i,2],extchg[i,3]))
out.append('$end')
for SectName, SectData in self.qctemplate:
for SectName, SectData, SectInsert in self.qctemplate:
if 'jobtype' in self.qcrems[remidx] and self.qcrems[remidx]['jobtype'].lower() == 'fsm':
fsm = True
if len(selection) != 2:
Expand All @@ -4233,15 +4247,19 @@ def write_qcin(self, **kwargs):
if SectName == 'molecule':
if molecule_printed == False:
molecule_printed = True
section_ln = 0
if read:
out.append("read")
section_ln = insert_section(out, section_ln, SectInsert)
elif self.na > 0:
out.append("%i %i" % (self.charge, self.mult))
section_ln = insert_section(out, section_ln, SectInsert)
an = 0
for e, x in zip(self.elem, self.xyzs[I]):
pre = '@' if ('qm_ghost' in self.Data and self.Data['qm_ghost'][an]) else ''
suf = self.Data['qcsuf'][an] if 'qcsuf' in self.Data else ''
out.append(pre + format_xyz_coord(e, x) + suf)
section_ln = insert_section(out, section_ln, SectInsert)
an += 1
if fsm:
out.append("****")
Expand All @@ -4250,6 +4268,7 @@ def write_qcin(self, **kwargs):
pre = '@' if ('qm_ghost' in self.Data and self.Data['qm_ghost'][an]) else ''
suf = self.Data['qcsuf'][an] if 'qcsuf' in self.Data else ''
out.append(pre + format_xyz_coord(e, x) + suf)
section_ln = insert_section(out, section_ln, SectInsert)
an += 1
if SectName == 'rem':
for key, val in self.qcrems[remidx].items():
Expand Down
Loading