Skip to content

Commit

Permalink
add more options to csearching
Browse files Browse the repository at this point in the history
  • Loading branch information
corinwagen committed Sep 4, 2020
1 parent c73cf56 commit 075a39e
Show file tree
Hide file tree
Showing 2 changed files with 17 additions and 12 deletions.
6 changes: 3 additions & 3 deletions cctk/molecule.py
Original file line number Diff line number Diff line change
Expand Up @@ -1648,14 +1648,14 @@ def optimize(self, inplace=True, nprocs=1):
else:
return optimized

def csearch(self, nprocs=1):
def csearch(self, nprocs=1, constraints=[]):
"""
Optimize molecule at the GFN2-xtb level of theory.
Args:
inplace (Bool): whether or not to return a new molecule or simply modify ``self.geometry``
nprocs (int): number of processors to use
constraints (list): atoms numbers to freeze
"""
import cctk.optimize as opt
assert isinstance(nprocs, int), "nprocs must be int!"
return opt.csearch(self, nprocs=nprocs)
return opt.csearch(self, nprocs=nprocs, constraints=constraints)
23 changes: 14 additions & 9 deletions cctk/optimize.py
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ def run_xtb(molecule, nprocs=1, return_energy=False):
else:
return output_mol

def csearch(molecule, constrain_atoms=None, rotamers=True, nprocs=1):
def csearch(molecule, constraints=None, rotamers=True, nprocs=1, noncovalent=False):
"""
Run a conformational search on a molecule using ``crest``.
Expand All @@ -87,6 +87,7 @@ def csearch(molecule, constrain_atoms=None, rotamers=True, nprocs=1):
constraints (list): list of atom numbers to constrain
rotamers (bool): return all rotamers or group into distinct conformers
nprocs (int): number of processors to use
noncovalent (Bool): whether or not to use non-covalent settings
Returns:
cctk.ConformationalEnsemble
Expand All @@ -97,19 +98,23 @@ def csearch(molecule, constrain_atoms=None, rotamers=True, nprocs=1):

assert installed("crest"), "crest must be installed!"

nci = ""
if noncovalent:
nci = "-nci"

ensemble = None
try:
if 1:
tmpdir = "/n/home03/cwagen/sw/cctk/crest/"
if constrain_atoms is not None:
assert isinstance(constrain_atoms, list)
assert all(isinstance(n, int) for n in constrain_atoms)
command = f"crest --constrain {','.join(constrain_atoms)}"
with tempfile.TemporaryDirectory() as tmpdir:
cctk.XYZFile.write_molecule_to_file(f"{tmpdir}/xtb-in.xyz", molecule)

if constraints is not None:
assert isinstance(constraints, list)
assert all(isinstance(n, int) for n in constraints)
command = f"crest xtb-in.xyz --constrain {','.join([str(c) for c in constraints])}"
result = sp.run(command, stdout=sp.PIPE, stderr=sp.PIPE, cwd=tmpdir, shell=True)
result.check_returncode()

cctk.XYZFile.write_molecule_to_file(f"{tmpdir}/xtb-in.xyz", molecule)
command = f"crest xtb-in.xyz --chrg {molecule.charge} --uhf {molecule.multiplicity - 1} -T {nprocs}"
command = f"crest xtb-in.xyz --chrg {molecule.charge} --uhf {molecule.multiplicity - 1} -T {nprocs} {nci}"
result = sp.run(command, stdout=sp.PIPE, stderr=sp.PIPE, cwd=tmpdir, shell=True)
result.check_returncode()

Expand Down

0 comments on commit 075a39e

Please sign in to comment.