merge conflix

cctk/molecule.py

@@ -1648,16 +1648,19 @@ def optimize(self, inplace=True, nprocs=1):
         else:
             return optimized
 
-    def csearch(self, nprocs=1):
+    def csearch(self, nprocs=1, constraints=[], logfile=None, noncovalent=False):
         """
         Optimize molecule at the GFN2-xtb level of theory.
 
         Args:
             nprocs (int): number of processors to use
+            constraints (list): atoms numbers to freeze
+            noncovalent (Bool): whether or not to use non-covalent settings
+            logfile (str): file to write ongoing ``crest`` output to
 
         Returns
             ConformationalEnsemble
         """
         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, noncovalent=noncovalent, logfile=logfile)

cctk/optimize.py

@@ -3,7 +3,7 @@
 """
 
 import numpy as np
-import os, tempfile, shutil, re
+import os, tempfile, shutil, re, logging
 import cctk
 import subprocess as sp
 
@@ -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=False, nprocs=1, noncovalent=False, logfile=None):
     """
     Run a conformational search on a molecule using ``crest``.
 
@@ -87,30 +87,42 @@ 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
+        logfile (str): file to write ongoing ``crest`` output to
 
     Returns:
         cctk.ConformationalEnsemble
     """
-
     assert isinstance(molecule, cctk.Molecule), "need a valid molecule!"
     assert isinstance(nprocs, int)
+    assert isinstance(logfile, str)
 
     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}"
-            result = sp.run(command, stdout=sp.PIPE, stderr=sp.PIPE, cwd=tmpdir, shell=True)
+            command = f"crest xtb-in.xyz --chrg {molecule.charge} --uhf {molecule.multiplicity - 1} -T {nprocs} {nci}"
+
+            if logfile:
+                with open(logfile, "w") as f:
+                    result = sp.run(command, stdout=f, stderr=f, cwd=tmpdir, shell=True)
+            else:
+                result = sp.run(command, stdout=sp.PIPE, stderr=sp.PIPE, cwd=tmpdir, shell=True)
+
             result.check_returncode()
 
             if rotamers: