Merge pull request #2 from AdaptivePELE/master

Adaptive 1.7.1 Release

AdaptivePELE/__init__.py

@@ -1,2 +1,2 @@
-__version__ = "1.7"
+__version__ = "1.7.1"
 name = "AdaptivePELE"

AdaptivePELE/adaptiveSampling.py

@@ -12,6 +12,7 @@
 import argparse
 import numpy as np
 from builtins import range
+import AdaptivePELE
 from AdaptivePELE.constants import blockNames, constants
 from AdaptivePELE.atomset import atomset
 from AdaptivePELE.utilities import utilities
@@ -610,7 +611,7 @@ def buildNewClusteringAndWriteInitialStructuresInNewSimulation(debug, controlFil
 
 def getClusteringLigandInfo(clustering_block_json):
     """
-        Get the ligand information (resname, resnum and chain from the cluering
+        Get the ligand information (resname, resnum and chain from the clustering
         parameters block in the control file)
 
         :param clustering_block_json: JSON block with the clustering-related parameters
@@ -620,9 +621,10 @@ def getClusteringLigandInfo(clustering_block_json):
             molecule to use for the clutering
 
     """
-    resname = str(clustering_block_json.get(blockNames.ClusteringTypes.ligandResname, "")).upper()
-    resnum = int(clustering_block_json.get(blockNames.ClusteringTypes.ligandResnum, 0))
-    resChain = str(clustering_block_json.get(blockNames.ClusteringTypes.ligandChain, "")).upper()
+    paramsBlock = clustering_block_json[blockNames.ClusteringTypes.params]
+    resname = str(paramsBlock.get(blockNames.ClusteringTypes.ligandResname, "")).upper()
+    resnum = int(paramsBlock.get(blockNames.ClusteringTypes.ligandResnum, 0))
+    resChain = str(paramsBlock.get(blockNames.ClusteringTypes.ligandChain, "")).upper()
     return resname, resnum, resChain
 
 def main(jsonParams, clusteringHook=None):
@@ -633,6 +635,7 @@ def main(jsonParams, clusteringHook=None):
         :type jsonParams: str
     """
 
+    utilities.print_unbuffered("Running AdaptivePELE version %s from %s" % (AdaptivePELE.__version__, AdaptivePELE.__path__[0]))
     controlFileValidator.validate(jsonParams)
     generalParams, spawningBlock, simulationrunnerBlock, clusteringBlock = loadParams(jsonParams)
 
@@ -776,7 +779,7 @@ def main(jsonParams, clusteringHook=None):
 
             degeneracyOfRepresentatives = spawningCalculator.calculate(clustersList, simulationRunner.getWorkingProcessors(), i, outputPathConstants=outputPathConstants)
             spawningCalculator.log()
-            # this method only does works with MSM-based spwaning methods,
+            # this method only does works with MSM-based spawning methods,
             # creating a plot of the stationary distribution and the PMF, for
             # the rest of methods it does nothing
             spawningCalculator.createPlots(outputPathConstants, i, clusteringMethod)

AdaptivePELE/analysis/addRejectedSteps.py

@@ -7,7 +7,7 @@
 """
 from __future__ import absolute_import, division, print_function, unicode_literals
 import os
-import glob
+from AdaptivePELE.utilities import utilities
 
 # to be changed by user
 reportFilename = "report_"
@@ -21,7 +21,7 @@
     print("Writing reports in folder", folder)
     os.chdir(folder)
 
-    reports = glob.glob(reportFilename + "*")
+    reports = utilities.getReportList(reportFilename + "*")
     for report in reports:
         with open(report, "r") as f:
             newReportContent = []

AdaptivePELE/analysis/bestStructs.py

@@ -1,12 +1,13 @@
 #!/usr/bin/python2.7
 
 import os
+import re
+import sys
+import glob
 import errno
 import argparse
 import pandas as pd
-import glob
-import re
-import sys
+from AdaptivePELE.utilities import utilities
 from AdaptivePELE.analysis import splitTrajectory
 
 """
@@ -82,8 +83,8 @@ def main(criteria, path=DIR, n_structs=10, sort_order="min", out_freq=FREQ, outp
         f_out: Name of the n outpu
     """
     #Get reports files
-    reports = glob.glob(os.path.join(path, "*/*report*"))
-    reports = glob.glob(os.path.join(path, "*report*")) if not reports else reports
+    reports = utilities.getReportList(os.path.join(path, "*/*report*"))
+    reports = utilities.getReportList(os.path.join(path, "*report*")) if not reports else reports
     reports = filter_non_numerical_folders(reports, numfolders)
     try:
         reports[0]

AdaptivePELE/analysis/cluster.py

@@ -19,6 +19,7 @@
 import re
 import sys
 import warnings
+from AdaptivePELE.utilities import utilities
 import AdaptivePELE.analysis.splitTrajectory as st
 import AdaptivePELE.analysis.backtrackAdaptiveTrajectory as bk
 
@@ -191,8 +192,8 @@ def Filter(values, percentage, column_name, thresh=None):
 
 
 def find_reports(path, numfolders):
-    reports = glob.glob(os.path.join(path, "*/*report_*"))
-    reports = glob.glob(os.path.join(path, "*report_*")) if not reports else reports
+    reports = utilities.getReportList(os.path.join(path, "*/*report*"))
+    reports = utilities.getReportList(os.path.join(path, "*report*")) if not reports else reports
     reports = filter_non_numerical_folders(reports, numfolders)
     try:
         reports[0]

AdaptivePELE/analysis/numberOfClusters.py

@@ -225,6 +225,8 @@ def plotContactsHistogram(folder, templetizedClusteringSummaryFile):
     """
     allFolders = os.listdir(folder)
     lastEpoch = len([epoch for epoch in allFolders if epoch.isdigit() and os.path.isfile(templetizedClusteringSummaryFile % int(epoch))]) - 1
+    if lastEpoch == -1:
+        raise ValueError("No simulation data in "+os.path.abspath(folder))
     lastSummary = templetizedClusteringSummaryFile % lastEpoch
     contactsColumn = 3
     allContacts = np.loadtxt(lastSummary, usecols=(contactsColumn,), ndmin=1)

AdaptivePELE/analysis/plotAdaptive.py

@@ -1,6 +1,5 @@
 from __future__ import absolute_import, division, print_function, unicode_literals
 import os
-import glob
 import argparse
 import matplotlib.pyplot as plt
 from AdaptivePELE.utilities import utilities
@@ -14,10 +13,7 @@ def parseArguments():
     """
         Parse command line arguments
 
-        :returns: int, int, int, str, bool, bool, int, str, bool, str, str, str, str -- Number of steps per epoch,
-            column to plot in the X axis, column to plot in the Y axis, name of
-            the files containing the simulation data, whether to plot the data
-            as points, wether to plot the data as lines, column to use as color, range of trajectories to select, whether to color each trajectory differently, path where to store the plot, label of the x-axis, label of the y-axis, label for the colorbar
+        :returns: argparse.Namespace -- Namespace with the input parameters
     """
     desc = "Plot relevant information from a simulation's report files.\n"\
            "It MUST be run from the root epoch folder (i.e., where it can find the folders 0/, 1/, 2/, ... lastEpoch/"
@@ -40,8 +36,11 @@ def parseArguments():
     parser.add_argument("--simulation_path", type=str, default=".", help="Path to the simulation output, defult is the current directory")
     parser.add_argument("--skip_first_step", action="store_true", help="Avoid plotting the first step of each report")
     parser.add_argument("--skip_steps", default=None, type=int, help="Number of steps to skip")
-    args = parser.parse_args()
-    return args.steps, args.xcol, args.ycol, args.filename, args.points, args.lines, args.zcol, args.traj_range, args.traj_col, args.output_path, args.xlabel, args.ylabel, args.cblabel, args.figure_size, args.show_plots, args.simulation_path, args.skip_first_step, args.skip_steps
+    parser.add_argument("--x_left_value", default=None, type=float, help="Value for the left end of the x axis")
+    parser.add_argument("--x_right_value", default=None, type=float, help="Value for the right end of the x axis")
+    parser.add_argument("--y_top_value", default=None, type=float, help="Value for the top end of the x axis")
+    parser.add_argument("--y_bottom_value", default=None, type=float, help="Value for the bottom end of the x axis")
+    return parser.parse_args()
 
 
 def addLine(data_plot, traj_num, epoch, steps, opt_dict, artists):
@@ -80,7 +79,8 @@ def addLine(data_plot, traj_num, epoch, steps, opt_dict, artists):
 def createPlot(reportName, column1, column2, stepsPerRun, printWithLines,
                paletteModifier, trajs_range=None, label_x=None, label_y=None,
                label_colorbar=None, fig_size=(6, 6), simulation_path=".",
-               skip_first_step=False, skip_steps=None):
+               skip_first_step=False, skip_steps=None, y_top=None,
+               y_bottom=None, x_left=None, x_right=None):
     """
         Generate a string to be passed to gnuplot
 
@@ -110,6 +110,14 @@ def createPlot(reportName, column1, column2, stepsPerRun, printWithLines,
         :type skip_first_step: bool
         :param skip_steps: Number of steps to skip in the plot
         :type skip_steps: int
+        :param y_bottom: Bottom limit of the y axis
+        :type y_bottom: float
+        :param y_top: Top limit of the y axis
+        :type y_top: float
+        :param x_left: Left limit of the x axis
+        :type x_bottom: float
+        :param x_right: Right limit of the x axis
+        :type x_right: float
     """
     epochs = utilities.get_epoch_folders(simulation_path)
     numberOfEpochs = int(len(epochs))
@@ -123,7 +131,7 @@ def createPlot(reportName, column1, column2, stepsPerRun, printWithLines,
     artists = []
     trajectory_range = set()
     if trajs_range is not None:
-        start, end = map(int, traj_range.split(":"))
+        start, end = map(int, trajs_range.split(":"))
         trajectory_range = set(range(start, end+1))
     cmin = 1e10
     cmax = -1e10
@@ -132,12 +140,10 @@ def createPlot(reportName, column1, column2, stepsPerRun, printWithLines,
     min_report = 1e10
     for epoch in epochs:
         ep = int(epoch)
-        reports = glob.glob(os.path.join(simulation_path, epoch, reportName+"*"))
+        reports = utilities.getReportList(os.path.join(simulation_path, epoch, reportName+"*"))
         if not reports:
             raise ValueError("Could not find any reports with the given name!!")
         for report in reports:
-            if not utilities.isReport(report):
-                continue
             report_num = utilities.getReportNum(report)
             max_report = max(max_report, report_num)
             min_report = min(min_report, report_num)
@@ -187,6 +193,8 @@ def createPlot(reportName, column1, column2, stepsPerRun, printWithLines,
             cbar.set_label("Epoch")
         if label_colorbar is not None:
             cbar.set_label(label_colorbar)
+    ax.set_ylim(bottom=y_bottom, top=y_top)
+    ax.set_xlim(left=x_left, right=x_right)
 
     annot = ax.annotate("", xy=(0, 0), xytext=(20, 20),
                         textcoords="offset points",
@@ -251,7 +259,7 @@ def hover(event):
 def generatePlot(stepsPerRun, xcol, ycol, reportName, kindOfPrint, paletteModifier,
                  trajs_range, path_to_save, xlabel, ylabel, cblabel, fig_size=(6, 6),
                  show_plot=True, simulation_path=".", skip_first_step=False,
-                 skip_steps=None):
+                 skip_steps=None, y_top=None, y_bottom=None, x_left=None, x_right=None):
     """
         Generate a template string to use with gnuplot
 
@@ -287,6 +295,16 @@ def generatePlot(stepsPerRun, xcol, ycol, reportName, kindOfPrint, paletteModifi
         :type skip_first_step: bool
         :param skip_first_step: Whether to avoid plotting the first point in each report
         :type skip_first_step: bool
+        :param skip_steps: Number of steps to skip in the plot
+        :type skip_steps: int
+        :param y_bottom: Bottom limit of the y axis
+        :type y_bottom: float
+        :param y_top: Top limit of the y axis
+        :type y_top: float
+        :param x_left: Left limit of the x axis
+        :type x_bottom: float
+        :param x_right: Right limit of the x axis
+        :type x_right: float
 
         :returns: str -- String to plot using gnuplot
     """
@@ -297,7 +315,8 @@ def generatePlot(stepsPerRun, xcol, ycol, reportName, kindOfPrint, paletteModifi
     createPlot(reportName, xcol, ycol, stepsPerRun, printWithLines, paletteModifier,
                trajs_range=trajs_range, label_x=xlabel, label_y=ylabel, label_colorbar=cblabel,
                fig_size=fig_size, simulation_path=simulation_path,
-               skip_first_step=skip_first_step, skip_steps=skip_steps)
+               skip_first_step=skip_first_step, skip_steps=skip_steps,
+               y_top=y_top, y_bottom=y_bottom, x_left=x_left, x_right=x_right)
     if path_to_save is not None:
         folder, _ = os.path.split(path_to_save)
         if folder:
@@ -307,8 +326,15 @@ def generatePlot(stepsPerRun, xcol, ycol, reportName, kindOfPrint, paletteModifi
         plt.show()
 
 if __name__ == "__main__":
-    steps_Run, Xcol, Ycol, filename, be, rmsd, colModifier, traj_range, color_traj, output_path, xlab, ylab, cblab, figure_size, plots_show, path_simulation, should_skip_first_step, n_skip_steps = parseArguments()
+    args = parseArguments()
+    be = args.points
+    rmsd = args.lines
+    colModifier = args.zcol
+    color_traj = args.traj_col
+    figure_size = args.figure_size
     figure_size = tuple(map(int, figure_size.split("x")))
+    Xcol = args.xcol
+    Ycol = args.ycol
     Xcol -= 1
     Ycol -= 1
     if colModifier is not None:
@@ -321,7 +347,9 @@ def generatePlot(stepsPerRun, xcol, ycol, reportName, kindOfPrint, paletteModifi
     if color_traj:
         colModifier = -1
 
-    generatePlot(steps_Run, Xcol, Ycol, filename, kind_Print, colModifier,
-                 traj_range, output_path, xlab, ylab, cblab, figure_size,
-                 show_plot=plots_show, simulation_path=path_simulation,
-                 skip_first_step=should_skip_first_step, skip_steps=n_skip_steps)
+    generatePlot(args.steps, Xcol, Ycol, args.filename, kind_Print, colModifier,
+                 args.traj_range, args.output_path, args.xlabel, args.ylabel,
+                 args.cblabel, figure_size, show_plot=args.show_plots,
+                 simulation_path=args.simulation_path, skip_first_step=args.skip_first_step,
+                 skip_steps=args.skip_steps, y_top=args.y_top_value, y_bottom=args.y_bottom_value,
+                 x_left=args.x_left_value, x_right=args.x_right_value)

AdaptivePELE/analysis/selectOnPlot.py

@@ -159,11 +159,7 @@ def concat_reports_in_csv(adaptive_results_path, output_file_path, report_prefix
     for adaptive_epoch in range(0, 2000):
         folder = os.path.join(adaptive_results_path, str(adaptive_epoch))
         if os.path.exists(folder):
-            report_list = glob.glob("{}/*{}*".format(folder, report_prefix))
-            # filter names that identify proper PELE reports, not interstep
-            # reports
-            report_list = [x for x in report_list if adapt_tools.isReport(x)]
-            report_list = sorted(report_list, key=adapt_tools.getReportNum)
+            report_list = adapt_tools.getReportList("{}/*{}*".format(folder, report_prefix))
             for n, report in enumerate(report_list):
                 if skip_first:
                     pandas_df = pd.read_csv(report, sep="    ", engine="python", index_col=False, header=0, skiprows=[1])

AdaptivePELE/analysis/simulationToCsv.py

@@ -18,8 +18,8 @@ def retrieve_fields(report):
     return list(data)[2:]
 
 def gather_reports():
-    reports = glob.glob(os.path.join("*/*report*"))
-    reports = glob.glob(os.path.join("*report*")) if not reports else reports
+    reports = utilities.getReportList(os.path.join(path, "*/*report*"))
+    reports = utilities.getReportList(os.path.join(path, "*report*")) if not reports else reports
     reports = filter_non_numerical_folders(reports)
     try:
         reports[0]

AdaptivePELE/clustering/clustering.py

@@ -1988,11 +1988,11 @@ def getAllTrajectories(paths):
 
         :param paths: The path where to find the trajectories
         :type paths: str
-        :returns: list -- A list with the names of all th trajectories in paths
+        :returns: list -- A list with the names of all the trajectories in paths
     """
     files = []
     for path in paths:
-        files += glob.glob(path)
+        files += utilities.getReportList(path)
     # sort the files obtained by glob by name, so that the results will be the
     # same on all computers
     return sorted(files)

AdaptivePELE/constants/blockNames.py

@@ -150,6 +150,9 @@ class SimulationParams:
     equilibrationSelect = "equilibrationSelect"
     equilibrationCluster = "equilibrationCluster"
     numberEquilibrationStructures = "numberEquilibrationStructures"
+    equilibrationRotationRange = "equilibrationRotationRange"
+    equilibrationTranslationRange = "equilibrationTranslationRange"
+    equilibrationBoxRadius = "equilibrationBoxRadius"
     boxCenter = "boxCenter"
     boxRadius = "boxRadius"
     runEquilibration = "runEquilibration"

AdaptivePELE/docs/Changelog.rst

@@ -10,6 +10,24 @@ The format is based on `Keep a Changelog <http://keepachangelog.com/en/1.0.0/>`_
 .. XXX - Unreleased
 .. ----------------
 
+1.7.1 - Unreleased
+------------------
+
+New features:
+.............
+
+    - Add equilibrationRotationRange, equilibrationTranslationRange and
+      equilibrationBoxRadius parameters to tune equilibration
+      performance
+    - Add support to correctly run simulations with PELE interstep reports and trajectories
+
+Bug fixes:
+..........
+
+    - Fix bug reading the ligand information in the adaptiveSampling main,
+      which caused to set the equilibration box in the center of mass of the
+      whole system
+
 1.7 - 2021-02-11
 ----------------
 

AdaptivePELE/docs/UserManual.rst

@@ -13,7 +13,7 @@ To install from PyPI simply run::
 
 To install from Conda run::
 
-     conda install -c nostrumbiodiscovery adaptive_pele 
+     conda install -c conda-forge -c nostrumbiodiscovery adaptive_pele 
 
 To install from source, you need to install and compile cython files in the base folder with::
 
@@ -190,6 +190,12 @@ Optionally, you can also use the following parameters:
 * **numberEquilibrationStructures** (*int*, default=10): Number of clusters to
   obtain from the *equilibrationCluster* structure selection (see
   **equilibrationMode** for more details)
+* **equilibrationBoxRadius** (*float*, default=2.0): Value of the simulation box 
+  radius for the equilibration (in angstroms)
+* **equilibrationRotationRange** (*float*, default=0.05): Value of the rotation 
+  magnitude in the equilibration simulation
+* **equilibrationTranslationRange** (*float*, default=0.5): Value of the
+  translation magnitude in the equilibration simulation
 * **useSrun** (*bool*, default=False): Whether to use srun to launch the PELE
   simulation instead of mpirun. Using srun allows a finer control over the
   resources used and might be helpful to deal with different cluster