replica exchange simulation submission + minor reeds tweaks

pygromos/files/blocks/imd_blocks.py

@@ -2,6 +2,7 @@
 
 from pygromos.files.blocks._general_blocks import TITLE as generic_TITLE
 from pygromos.files.blocks._general_blocks import _generic_gromos_block
+from pygromos.utils.utils import str2bool
 from pygromos.utils.typing import List, Dict, Union, Number
 
 # forward declarations
@@ -322,31 +323,31 @@ def __init__(
     ):
         super().__init__(used=True, content=content)
         if content is None:
-            self.RETL = bool(RETL)
+            self.RETL = str2bool(RETL)
 
             self.NRET = int(NRET)
             self.RET = RET
 
-            self.LRESCALE = bool(LRESCALE)
+            self.LRESCALE = str2bool(LRESCALE)
 
             self.NRELAM = int(NRELAM)
             self.RELAM = RELAM
             self.RETS = RETS
 
             self.NRETRIAL = int(NRETRIAL)
             self.NREQUIL = int(NREQUIL)
-            self.CONT = bool(CONT)
+            self.CONT = str2bool(CONT)
 
     def read_content_from_str(self, content: List[str]):
         try:
-            setattr(self, "RETL", int(content[1].split()[0]))
+            setattr(self, "RETL", str2bool(content[1].split()[0]))
             setattr(self, "NRET", int(content[3].split()[0]))
             T_values = list(map(float, content[5].split()))
             if len(T_values) == self.NRET:
                 setattr(self, "RET", T_values)
             else:
                 raise IOError("REPLICA: NRET was not equal to the number of Temperatures (RET) in IMD!")
-            setattr(self, "LRESCALE", int(content[7].split()[0]))
+            setattr(self, "LRESCALE", str2bool(content[7].split()[0]))
             setattr(self, "NRELAM", int(content[9].split()[0]))
             lambda_val = list(map(float, content[11].split()))
             if len(lambda_val) == self.NRELAM:
@@ -363,7 +364,7 @@ def read_content_from_str(self, content: List[str]):
             raise IOError("Could not parse block from str - " + __class__.__name__ + "\n" + str(err.args))
 
 
-class NEW_REPLICA_EDS(_generic_imd_block):
+class REPLICA_EDS(_generic_imd_block):
     """REPLICA_EDS Block
 
         This block is controlling the REPLICA_EDS settings in  gromos and is basically a mixture of EDS and RE block. (Don't use them when using this block!)
@@ -406,8 +407,8 @@ class NEW_REPLICA_EDS(_generic_imd_block):
 
     NRETRIAL: int
     NREQUIL: int
-    EDS_STAT_OUT: int
     CONT: bool
+    EDS_STAT_OUT: int
     PERIODIC: int
 
     _order = [
@@ -437,7 +438,7 @@ def __init__(
     ):
         super().__init__(used=True, content=content)
         if content is None:
-            self.REEDS = bool(REEDS)
+            self.REEDS = str2bool(REEDS)
 
             self.NRES = int(NRES)
             self.NEOFF = int(NEOFF)
@@ -448,119 +449,16 @@ def __init__(
 
             self.NRETRIAL = int(NRETRIAL)
             self.NREQUIL = int(NREQUIL)
-            self.CONT = bool(CONT)
+            self.CONT = str2bool(CONT)
             self.EDS_STAT_OUT = int(EDS_STAT_OUT)
             self.PERIODIC = int(PERIODIC)
 
     def read_content_from_str(self, content: List[str]):
         try:
-            setattr(self, "REEDS", int(content[1].split()[0]))
-            setattr(self, "NRES", int(content[3].split()[0]))
-            setattr(self, "NEOFF", int(content[3].split()[1]))
-            setattr(self, "NUMSTATES", int(content[3].split()[2]))
-            s_values = list(map(float, content[5].split()))
-            if len(s_values) == self.NRES:
-                setattr(self, "RES", s_values)
-            else:
-                raise IOError("REPLICA_EDS: NRES was not equal to the number of s-values in IMD!")
-            EIR = []
-            for ind in range(7, 7 + self.NUMSTATES):
-                EIR_line = list(map(float, content[ind].split()))
-                if len(EIR_line) != self.NRES:
-                    raise IOError("REPLICA_EDS: NRES was not equal to the number of EIRs given in IMD!")
-                EIR.append(EIR_line)
-            setattr(self, "EIR", EIR)
-            [setattr(self, key, int(value)) for key, value in zip(self._order[0][-1], content[-1].split())]
-
-        except Exception as err:
-            raise IOError("Could not parse block from str - " + __class__.__name__ + "\n" + str(err.args))
-
-
-class REPLICA_EDS(_generic_imd_block):
-    name: str = "REPLICA_EDS"
-
-    REEDS: bool
-
-    NRES: int
-    NUMSTATES: int
-
-    RES: List[float]
-    EIR: List[float]
-
-    NRETRIAL: int
-    NREQUIL: int
-    EDS_STAT_OUT: int
-    CONT: bool
-
-    _order = [
-        [
-            ["REEDS"],
-            ["NRES", "NUMSTATES"],
-            ["RES(1 ... NRES)"],
-            ["EIR(NUMSTATES x NRES)"],
-            ["NRETRIAL", "NREQUIL", "CONT", "EDS_STAT_OUT"],
-        ]
-    ]
-
-    def __init__(
-        self,
-        REEDS: bool = True,
-        NRES: int = 0,
-        NUMSTATES: int = 0,
-        RES: List[float] = [0],
-        EIR: List[List[float]] = [[0]],
-        NRETRIAL: int = 0,
-        NREQUIL: int = 0,
-        EDS_STAT_OUT: int = 0,
-        CONT: bool = True,
-        content: List[str] = None,
-    ):
-        """REPLICA_EDS Block
-
-            This block is controlling the REPLICA_EDS settings in  gromos and is basically a mixture of EDS and RE block. (Don't use them when using this block!)
-
-        Attributes
-        ----------
-        REEDS:  bool
-            Shall REEDS be activated?
-        NRES:   int
-            Number of s-Values
-        NUMSTATES:  int
-            Number of EDS-states
-
-        RES:    List[float]
-            s_values for all replicas
-        EIR:    List[List[float]]
-            energy offsets for all replicas and all states  List[List[float]] = REPLICA[EDS_STATE[EIR]]
-        NERTRIAL: int
-            How many replica exchanges trials should be executed? (NRETRIAL*STEP.NSTLIM == total simulation time)
-        NREQUIL: int
-            How many equilibration runs shall be exectured? (NREQUIL*STEP.NSTLIM == total simulation time)
-        EDS_STAT_OUT: int
-            Shall the replica exchange information be outputted? (__future__ frequency of output.)
-        CONT: bool
-            Is this a continuation run?
-        """
-        super().__init__(used=True, content=content)
-        if content is None:
-            self.REEDS = REEDS
-
-            self.NRES = NRES
-            self.NUMSTATES = NUMSTATES
-
-            self.RES = RES
-            self.EIR = EIR
-
-            self.NRETRIAL = NRETRIAL
-            self.NREQUIL = NREQUIL
-            self.CONT = CONT
-            self.EDS_STAT_OUT = EDS_STAT_OUT
-
-    def read_content_from_str(self, content: List[str]):
-        try:
-            setattr(self, "REEDS", int(content[1].split()[0]))
+            setattr(self, "REEDS", str2bool(content[1].split()[0]))
             setattr(self, "NRES", int(content[3].split()[0]))
             setattr(self, "NUMSTATES", int(content[3].split()[1]))
+            setattr(self, "NEOFF", int(content[3].split()[2]))
             s_values = list(map(float, content[5].split()))
             if len(s_values) == self.NRES:
                 setattr(self, "RES", s_values)

pygromos/files/simulation_parameters/imd.py

@@ -58,7 +58,6 @@ class Imd(_general_gromos_file._general_gromos_file):
     INITIALISE: blocks.INITIALISE
 
     REPLICA_EDS: blocks.REPLICA_EDS
-    NEW_REPLICA_EDS: blocks.NEW_REPLICA_EDS
 
     REPLICA: blocks.REPLICA
 
@@ -109,7 +108,7 @@ def edit_EDS(
             self.EDS.EIR = EIR
 
         else:
-            print("Setting new EDS_block")
+            # print("Setting new EDS_block")
             if type(EIR) == float or type(EIR) == str or type(EIR) == int:
                 EIR = [float(EIR) for x in range(NUMSTATES)]
 
@@ -133,7 +132,7 @@ def edit_REEDS(
 
         # specific relations are rescued here
         reeds_block = self.REPLICA_EDS
-        print(type(reeds_block))
+        # print(type(reeds_block))
 
         if isinstance(REEDS, bool):
             reeds_block.REEDS = REEDS

pygromos/simulations/hpc_queuing/job_scheduling/workers/simulation_workers/simulation_run_worker.py

@@ -132,7 +132,24 @@ def work(
     tmp_prefix = os.path.basename(out_dir)
     tmp_imd_path = out_dir + "/" + tmp_prefix + ".imd"
     imd_file = imd.Imd(in_imd_path)
-    cnf_file = cnf.Cnf(in_cnf_path)
+
+    if hasattr(imd_file, "REPLICA") and imd_file.REPLICA is not None:
+        is_repex_run = True
+        num_replicas = int(imd_file.REPLICA.NRELAM * imd_file.REPLICA.NRET)
+    elif hasattr(imd_file, "REPLICA_EDS") and imd_file.REPLICA_EDS is not None:
+        is_repex_run = True
+        num_replicas = int(imd_file.REPLICA_EDS.NRES)
+    else:
+        is_repex_run = False
+
+    # Prepare CNF file(s):
+
+    if is_repex_run:
+        tmp_path = "/".join(os.path.abspath(in_cnf_path).split("/")[:-1])
+        in_coord_paths = sorted(glob.glob(tmp_path + "/*.cnf"))
+        cnf_file = cnf.Cnf(in_coord_paths[0])  # used to make sure imd/cnf is in sync
+    else:
+        cnf_file = cnf.Cnf(in_cnf_path)
 
     # check init_block - if specified!
     # What kind of simulation
@@ -178,6 +195,12 @@ def work(
         if is_stochastic_dynamics_sim or is_vacuum:
             imd_file.INITIALISE.NTISTI = 1
 
+    # For RE-EDS / REPEX simulations, also change the CONT option to 1
+    if hasattr(imd_file, "REPLICA") and runID > 1:
+        imd_file.REPLICA.CONT = 1
+    if hasattr(imd_file, "REPLICA_EDS") and runID > 1:
+        imd_file.REPLICA_EDS.CONT = 1
+
     # adjust sim time if continuation:
     if runID > 1:
         imd_file.STEP.T = imd_file.STEP.T + (imd_file.STEP.NSTLIM * imd_file.STEP.DT) * (runID - 1)
@@ -195,35 +218,66 @@ def work(
         print(spacer + "\n start MD " + str(os.path.basename(tmp_imd_path)) + "\n")
 
         # TODO: This is a stupid workaround as Euler tends to place nans in the euler angles, that should not be there!
-        if hasattr(cnf_file, "GENBOX") and any([math.isnan(x) for x in cnf_file.GENBOX.euler]):
+        if is_repex_run:  # do the workaround for each cnf one by one
+            for tmp_cnf_path in in_coord_paths:
+                tmp_cnf = cnf.Cnf(tmp_cnf_path)
+                if hasattr(tmp_cnf, "GENBOX") and any([math.isnan(x) for x in tmp_cnf.GENBOX.euler]):
+                    tmp_cnf.GENBOX.euler = [0.0, 0.0, 0.0]
+                    tmp_cnf.write(tmp_cnf_path)
+        elif hasattr(cnf_file, "GENBOX") and any([math.isnan(x) for x in cnf_file.GENBOX.euler]):
             cnf_file.GENBOX.euler = [0.0, 0.0, 0.0]
             cnf_file.write(in_cnf_path)
+
+        # Start the execution of the gromosXX binary
         try:
-            omd_file_path = gromosXX.md_run(
-                in_topo_path=in_top_path,
-                in_coord_path=in_cnf_path,
-                in_imd_path=tmp_imd_path,
-                in_pert_topo_path=in_perttopo_path,
-                in_disres_path=in_disres_path,
-                in_posresspec_path=in_posres_path,
-                in_refpos_path=in_refpos_path,
-                in_qmmm_path=in_qmmm_path,
-                nmpi=nmpi,
-                nomp=nomp,
-                out_prefix=tmp_prefix,
-                out_tre=out_tre,
-                out_trc=out_trc,
-                out_trg=out_trg,
-                out_trs=out_trs,
-                out_trf=out_trf,
-                out_trv=out_trv,
-                verbose=True,
-            )
-
-            print("Waiting to find: ", omd_file_path.replace(".omd", ".cnf"))
-            bash.wait_for_fileSystem(omd_file_path.replace(".omd", ".cnf"))
-
-            md_failed = False
+            if is_repex_run:
+                omd_file_path = gromosXX.repex_run(
+                    in_topo_path=in_top_path,
+                    in_coord_path=in_cnf_path,
+                    in_imd_path=tmp_imd_path,
+                    in_pert_topo_path=in_perttopo_path,
+                    in_disres_path=in_disres_path,
+                    in_posresspec_path=in_posres_path,
+                    in_refpos_path=in_refpos_path,
+                    nmpi=nmpi,
+                    nomp=nomp,
+                    out_prefix=tmp_prefix,
+                    out_tre=out_tre,
+                    out_trc=out_trc,
+                    out_trg=out_trg,
+                    out_trs=out_trs,
+                    out_trf=out_trf,
+                    out_trv=out_trv,
+                    verbose=True,
+                )
+                md_failed = False
+
+            else:
+                omd_file_path = gromosXX.md_run(
+                    in_topo_path=in_top_path,
+                    in_coord_path=in_cnf_path,
+                    in_imd_path=tmp_imd_path,
+                    in_pert_topo_path=in_perttopo_path,
+                    in_disres_path=in_disres_path,
+                    in_posresspec_path=in_posres_path,
+                    in_refpos_path=in_refpos_path,
+                    in_qmmm_path=in_qmmm_path,
+                    nmpi=nmpi,
+                    nomp=nomp,
+                    out_prefix=tmp_prefix,
+                    out_tre=out_tre,
+                    out_trc=out_trc,
+                    out_trg=out_trg,
+                    out_trs=out_trs,
+                    out_trf=out_trf,
+                    out_trv=out_trv,
+                    verbose=True,
+                )
+
+                print("Waiting to find: ", omd_file_path.replace(".omd", ".cnf"))
+                bash.wait_for_fileSystem(omd_file_path.replace(".omd", ".cnf"))
+
+                md_failed = False
         except Exception as err:
             print("Failed! process returned: \n Err: \n" + "\n".join(err.args))
             md_failed = True
@@ -241,14 +295,30 @@ def work(
                 bash.move_file(work_dir + "/*", out_dir)
             else:
                 for host in hosts:
-                    command = "ssh " + host + '  "mv ${TMPDIR}/* ' + out_dir + '"'
+                    if host == os.environ["HOSTNAME"]:
+                        command = "mv ${TMPDIR}/* " + out_dir
+                    else:
+                        command = "ssh " + host + '  "mv ${TMPDIR}/* ' + out_dir + '"'
                     os.system(command)
             os.system("remote_tmpdir delete")  # Works for both multi or single node
 
         # Note: If job is multi-node, it is simpler to zip things in out_dir after copying back
         if multi_node and zip_trajectories:
             zip_files.do(in_simulation_dir=out_dir, n_processes=n_cpu_zip)
 
+        # Check for the output files:
+        if is_repex_run or multi_node:
+            try:
+                print("Ensuring we have all output cnfs: ", omd_file_path.replace(".omd", "*.cnf"))
+                out_cnf_paths = [omd_file_path.replace(".omd", "_" + str(n + 1) + ".cnf") for n in range(num_replicas)]
+                print(out_cnf_paths)
+                for ocp in out_cnf_paths:
+                    bash.wait_for_fileSystem(out_dir + "/" + ocp)
+            except Exception as err:
+                print("Failed! process returned: \n Err: \n" + "\n".join(err.args))
+                print("Missing CNF files.")
+                md_failed = True
+
     except Exception as err:
         print("\nFailed during simulations: ", file=sys.stderr)
         print(type(err), file=sys.stderr)

pygromos/utils/utils.py

@@ -31,6 +31,8 @@ def _cartesian_distance(x1: float, x2: float, y1: float, y2: float, z1: float, z
 def str2bool(v):
     if isinstance(v, bool):
         return v
+    elif isinstance(v, int) or isinstance(v, float):
+        return bool(v)
     if v.lower() in ("yes", "true", "t", "y", "1"):
         return True
     elif v.lower() in ("no", "false", "f", "n", "0"):