From e45c7b5ef17e4a14f8301b489d55c01481811d0c Mon Sep 17 00:00:00 2001 From: Maciej Bartkowiak Date: Mon, 18 Mar 2024 11:39:59 +0000 Subject: [PATCH 1/3] Correct unit cell detection in ASE converter --- MDANSE/Src/MDANSE/Framework/Converters/ASE.py | 12 +++++++++--- .../Framework/InputData/HDFTrajectoryInputData.py | 6 +----- 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py index 04534c9278..9719369f83 100644 --- a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py +++ b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py @@ -91,6 +91,7 @@ def initialize(self): """ Initialize the job. """ + self._fractionalCoordinates = None self._atomicAliases = self.configuration["atom_aliases"]["value"] # The number of steps of the analysis. @@ -101,6 +102,9 @@ def initialize(self): ).toval("ps") self.parse_first_step(self._atomicAliases) + print( + f"Fractional coords after parse_first_step: {self._fractionalCoordinates}" + ) self._start = 0 if self.numberOfSteps < 1: @@ -137,6 +141,7 @@ def run_step(self, index): time = self._timeaxis[index] unitCell = frame.cell.array + print(f"Unit cell from frame: {unitCell}") unitCell *= measure(1.0, "ang").toval("nm") unitCell = UnitCell(unitCell) @@ -145,10 +150,9 @@ def run_step(self, index): coords *= measure(1.0, "ang").toval("nm") if self._fractionalCoordinates: - conf = PeriodicBoxConfiguration( + realConf = PeriodicBoxConfiguration( self._trajectory.chemical_system, coords, unitCell ) - realConf = conf.to_real_configuration() else: realConf = PeriodicRealConfiguration( self._trajectory.chemical_system, coords, unitCell @@ -206,7 +210,9 @@ def parse_first_step(self, mapping): self._timeaxis = self._timestep * np.arange(self._total_number_of_steps) - self._fractionalCoordinates = np.all(first_frame.get_pbc()) + if self._fractionalCoordinates is None: + self._fractionalCoordinates = np.all(first_frame.get_pbc()) + print(f"PBC in first frame = {first_frame.get_pbc()}") g = Graph() diff --git a/MDANSE/Src/MDANSE/Framework/InputData/HDFTrajectoryInputData.py b/MDANSE/Src/MDANSE/Framework/InputData/HDFTrajectoryInputData.py index 20a61f9bc0..ca30d4de42 100644 --- a/MDANSE/Src/MDANSE/Framework/InputData/HDFTrajectoryInputData.py +++ b/MDANSE/Src/MDANSE/Framework/InputData/HDFTrajectoryInputData.py @@ -43,11 +43,7 @@ def info(self): val.append("Number of steps:") val.append("%s\n" % len(self._data)) val.append("Configuration:") - val.append( - "\tIs periodic: {}\n".format( - "unit_cell" in self._data.file["/configuration"] - ) - ) + val.append("\tIs periodic: {}\n".format("unit_cell" in self._data.file)) val.append("Variables:") for k, v in self._data.file["/configuration"].items(): val.append("\t- {}: {}".format(k, v.shape)) From 1adaea57f526ea4c3df34a1bc0536e5b68f4fafd Mon Sep 17 00:00:00 2001 From: Maciej Bartkowiak Date: Mon, 18 Mar 2024 13:07:02 +0000 Subject: [PATCH 2/3] Make MolecularViewer skip the unit cell iif it is not defined --- MDANSE/Src/MDANSE/Framework/Converters/ASE.py | 29 ++++---- .../MolecularViewer/MolecularViewer.py | 66 ++++++++++--------- .../MolecularViewer/readers/hdf5wrapper.py | 5 +- 3 files changed, 55 insertions(+), 45 deletions(-) diff --git a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py index 9719369f83..b166e94442 100644 --- a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py +++ b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py @@ -30,6 +30,7 @@ from MDANSE.MolecularDynamics.Configuration import ( PeriodicBoxConfiguration, PeriodicRealConfiguration, + RealConfiguration, ) from MDANSE.MolecularDynamics.Trajectory import TrajectoryWriter from MDANSE.MolecularDynamics.UnitCell import UnitCell @@ -91,7 +92,7 @@ def initialize(self): """ Initialize the job. """ - self._fractionalCoordinates = None + self._isPeriodic = None self._atomicAliases = self.configuration["atom_aliases"]["value"] # The number of steps of the analysis. @@ -102,9 +103,7 @@ def initialize(self): ).toval("ps") self.parse_first_step(self._atomicAliases) - print( - f"Fractional coords after parse_first_step: {self._fractionalCoordinates}" - ) + print(f"isPeriodic after parse_first_step: {self._isPeriodic}") self._start = 0 if self.numberOfSteps < 1: @@ -140,22 +139,24 @@ def run_step(self, index): frame = read(self.configuration["trajectory_file"]["value"], index=index) time = self._timeaxis[index] - unitCell = frame.cell.array - print(f"Unit cell from frame: {unitCell}") + if self._isPeriodic: + unitCell = frame.cell.array + print(f"Unit cell from frame: {unitCell}") - unitCell *= measure(1.0, "ang").toval("nm") - unitCell = UnitCell(unitCell) + unitCell *= measure(1.0, "ang").toval("nm") + unitCell = UnitCell(unitCell) coords = frame.get_positions() coords *= measure(1.0, "ang").toval("nm") - if self._fractionalCoordinates: - realConf = PeriodicBoxConfiguration( + if self._isPeriodic: + realConf = PeriodicRealConfiguration( self._trajectory.chemical_system, coords, unitCell ) else: - realConf = PeriodicRealConfiguration( - self._trajectory.chemical_system, coords, unitCell + realConf = RealConfiguration( + self._trajectory.chemical_system, + coords, ) self._trajectory.chemical_system.configuration = realConf @@ -210,8 +211,8 @@ def parse_first_step(self, mapping): self._timeaxis = self._timestep * np.arange(self._total_number_of_steps) - if self._fractionalCoordinates is None: - self._fractionalCoordinates = np.all(first_frame.get_pbc()) + if self._isPeriodic is None: + self._isPeriodic = np.all(first_frame.get_pbc()) print(f"PBC in first frame = {first_frame.get_pbc()}") g = Graph() diff --git a/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/MolecularViewer.py b/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/MolecularViewer.py index 888b691ef2..efc806b78f 100644 --- a/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/MolecularViewer.py +++ b/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/MolecularViewer.py @@ -646,36 +646,42 @@ def set_coordinates(self, frame: int, tolerance=0.04): if self._cell_visible: # update the unit cell uc = self._reader.read_pbc(self._current_frame) - a = uc.a_vector - b = uc.b_vector - c = uc.c_vector - uc_points = vtk.vtkPoints() - uc_points.SetNumberOfPoints(8) - for i, v in enumerate([[0, 0, 0], a, b, c, a + b, a + c, b + c, a + b + c]): - x, y, z = v - uc_points.SetPoint(i, x, y, z) - self._uc_polydata.SetPoints(uc_points) - - uc_lines = vtk.vtkCellArray() - for i, j in [ - (0, 1), - (0, 2), - (0, 3), - (1, 4), - (1, 5), - (4, 7), - (2, 4), - (2, 6), - (5, 7), - (3, 5), - (3, 6), - (6, 7), - ]: - line = vtk.vtkLine() - line.GetPointIds().SetId(0, i) - line.GetPointIds().SetId(1, j) - uc_lines.InsertNextCell(line) - self._uc_polydata.SetLines(uc_lines) + if uc is not None: + a = uc.a_vector + b = uc.b_vector + c = uc.c_vector + uc_points = vtk.vtkPoints() + uc_points.SetNumberOfPoints(8) + for i, v in enumerate( + [[0, 0, 0], a, b, c, a + b, a + c, b + c, a + b + c] + ): + x, y, z = v + uc_points.SetPoint(i, x, y, z) + self._uc_polydata.SetPoints(uc_points) + + uc_lines = vtk.vtkCellArray() + for i, j in [ + (0, 1), + (0, 2), + (0, 3), + (1, 4), + (1, 5), + (4, 7), + (2, 4), + (2, 6), + (5, 7), + (3, 5), + (3, 6), + (6, 7), + ]: + line = vtk.vtkLine() + line.GetPointIds().SetId(0, i) + line.GetPointIds().SetId(1, j) + uc_lines.InsertNextCell(line) + self._uc_polydata.SetLines(uc_lines) + else: + uc_lines = vtk.vtkCellArray() + self._uc_polydata.SetLines(uc_lines) else: uc_lines = vtk.vtkCellArray() self._uc_polydata.SetLines(uc_lines) diff --git a/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/readers/hdf5wrapper.py b/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/readers/hdf5wrapper.py index 7f4cb783c0..7b3d95cc72 100644 --- a/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/readers/hdf5wrapper.py +++ b/MDANSE_GUI/Src/MDANSE_GUI/MolecularViewer/readers/hdf5wrapper.py @@ -37,5 +37,8 @@ def read_frame(self, frame: int) -> "np.array": return np.array(coords) def read_pbc(self, frame: int) -> "np.array": - unit_cell = self._chemical_system.configuration.unit_cell + try: + unit_cell = self._chemical_system.configuration.unit_cell + except AttributeError: + unit_cell = None return unit_cell From 10754afb4a73a3cda65162677895cef7de65b029 Mon Sep 17 00:00:00 2001 From: Maciej Bartkowiak Date: Mon, 18 Mar 2024 13:28:54 +0000 Subject: [PATCH 3/3] Improve the performance of ASE converter --- MDANSE/Src/MDANSE/Framework/Converters/ASE.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py index b166e94442..cd5c031e37 100644 --- a/MDANSE/Src/MDANSE/Framework/Converters/ASE.py +++ b/MDANSE/Src/MDANSE/Framework/Converters/ASE.py @@ -136,6 +136,9 @@ def run_step(self, index): try: frame = self._input[index] except TypeError: + frame = next(self._input) + else: + print("ASE using the slower way") frame = read(self.configuration["trajectory_file"]["value"], index=index) time = self._timeaxis[index] @@ -195,15 +198,15 @@ def parse_first_step(self, mapping): try: self._input = ASETrajectory(self.configuration["trajectory_file"]["value"]) except: - self._input = iread( - self.configuration["trajectory_file"]["value"], index="[:]" - ) first_frame = read(self.configuration["trajectory_file"]["value"], index=0) last_iterator = 0 generator = iread(self.configuration["trajectory_file"]["value"]) for _ in generator: last_iterator += 1 generator.close() + self._input = iread( + self.configuration["trajectory_file"]["value"] # , index="[:]" + ) self._total_number_of_steps = last_iterator else: first_frame = self._input[0]