From e9e830a966a58304123f0a1d135db360c02caf58 Mon Sep 17 00:00:00 2001 From: jrudz Date: Wed, 4 Dec 2024 21:54:28 +0100 Subject: [PATCH 1/3] initial investigation into harmonizing particles and atoms --- .../schema_packages/atoms_state.py | 16 ++++- .../schema_packages/general.py | 58 +++++++++++-------- .../schema_packages/model_system.py | 16 ++--- .../schema_packages/particles_state.py | 10 +++- 4 files changed, 67 insertions(+), 33 deletions(-) diff --git a/src/nomad_simulations/schema_packages/atoms_state.py b/src/nomad_simulations/schema_packages/atoms_state.py index 32fbdd11..ed2c6526 100644 --- a/src/nomad_simulations/schema_packages/atoms_state.py +++ b/src/nomad_simulations/schema_packages/atoms_state.py @@ -552,7 +552,17 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: ) -class AtomsState(Entity): +class State(Entity): + """ + A base section to define the state information of the system. + """ + + def __init__(self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs): + super().__init__(m_def, m_context, **kwargs) + self.labels = None + + +class AtomsState(State): """ A base section to define each atom state information. """ @@ -641,3 +651,7 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: self.chemical_symbol = self.resolve_chemical_symbol(logger=logger) if self.atomic_number is None: self.atomic_number = self.resolve_atomic_number(logger=logger) + + # Set the labels + if self.chemical_symbol is not None: + self.labels = self.chemical_symbol diff --git a/src/nomad_simulations/schema_packages/general.py b/src/nomad_simulations/schema_packages/general.py index 59777858..496cebad 100644 --- a/src/nomad_simulations/schema_packages/general.py +++ b/src/nomad_simulations/schema_packages/general.py @@ -1,10 +1,17 @@ -from typing import TYPE_CHECKING +#! TODO: Why is TYPE_CHECKING False? +from typing import TYPE_CHECKING, List, Iterable, Union -if TYPE_CHECKING: - from collections.abc import Callable - - from nomad.datamodel.datamodel import EntryArchive - from structlog.stdlib import BoundLogger +if not TYPE_CHECKING: + from nomad.datamodel.datamodel import ( + EntryArchive, + ) + from nomad.metainfo import ( + Context, + Section, + ) + from structlog.stdlib import ( + BoundLogger, + ) import numpy as np from nomad.config import config @@ -227,7 +234,7 @@ def resolve_composition_formula(self, system_parent: ModelSystem) -> None: """ def set_composition_formula( - system: ModelSystem, subsystems: list[ModelSystem], atom_labels: list[str] + system: ModelSystem, subsystems: list[ModelSystem], labels: list[str] ) -> None: """Determine the composition formula for `system` based on its `subsystems`. If `system` has no children, the atom_labels are used to determine the formula. @@ -243,8 +250,8 @@ def set_composition_formula( system.atom_indices if system.atom_indices is not None else [] ) subsystem_labels = ( - [np.array(atom_labels)[atom_indices]] - if atom_labels + [np.array(labels)[atom_indices]] + if labels else ['Unknown' for atom in range(len(atom_indices))] ) else: @@ -259,7 +266,7 @@ def set_composition_formula( children_names=subsystem_labels ) - def get_composition_recurs(system: ModelSystem, atom_labels: list[str]) -> None: + def get_composition_recurs(system: ModelSystem, labels: list[str]) -> None: """Traverse the system hierarchy downward and set the branch composition for all (sub)systems at each level. @@ -269,23 +276,28 @@ def get_composition_recurs(system: ModelSystem, atom_labels: list[str]) -> None: to the atom indices stored in system. """ subsystems = system.model_system - set_composition_formula( - system=system, subsystems=subsystems, atom_labels=atom_labels - ) + set_composition_formula(system=system, subsystems=subsystems, labels=labels) if subsystems: for subsystem in subsystems: - get_composition_recurs(system=subsystem, atom_labels=atom_labels) + get_composition_recurs(system=subsystem, labels=labels) # ! CG: system_parent.cell[0].particles_state instead of atoms_state! - atoms_state = ( - system_parent.cell[0].atoms_state if system_parent.cell is not None else [] - ) - atom_labels = ( - [atom.chemical_symbol for atom in atoms_state] - if atoms_state is not None - else [] - ) - get_composition_recurs(system=system_parent, atom_labels=atom_labels) + labels = [] + if system_parent.cell is not None: + if system_parent.cell[0].name == 'AtomicCell': + labels = ( + [atom.labels for atom in system_parent.cell[0].atoms_state] + if system_parent.cell[0].atoms_state is not None + else [] + ) + elif system_parent.cell[0].name == 'ParticleCell': + labels = ( + [atom.labels for atom in system_parent.cell[0].particles_state] + if system_parent.cell[0].particles_state is not None + else [] + ) + + get_composition_recurs(system=system_parent, labels=labels) def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: super(Schema, self).normalize(archive, logger) diff --git a/src/nomad_simulations/schema_packages/model_system.py b/src/nomad_simulations/schema_packages/model_system.py index 21cd6ad4..d0cb0042 100644 --- a/src/nomad_simulations/schema_packages/model_system.py +++ b/src/nomad_simulations/schema_packages/model_system.py @@ -1372,10 +1372,12 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: #! ChemicalFormula calls `ase_atoms = atomic_cell.to_ase_atoms(logger=logger)` and `ase_atoms.get_chemical_formula()` # Creating and normalizing ChemicalFormula section - # TODO add support for fractional formulas (possibly add `AtomicCell.concentrations` for each species) - sec_chemical_formula = self.m_create(ChemicalFormula) - sec_chemical_formula.normalize(archive, logger) - if sec_chemical_formula.m_cache: - self.elemental_composition = sec_chemical_formula.m_cache.get( - 'elemental_composition', [] - ) + if any(cell.name == 'AtomicCell' for cell in self.cell): + # TODO: get_sibling_section() may need to be updated to more specifically search for AtomicCell in ChemicalFormula and Symmetry, in cases where multiple different cells are present + # TODO add support for fractional formulas (possibly add `AtomicCell.concentrations` for each species) + sec_chemical_formula = self.m_create(ChemicalFormula) + sec_chemical_formula.normalize(archive, logger) + if sec_chemical_formula.m_cache: + self.elemental_composition = sec_chemical_formula.m_cache.get( + 'elemental_composition', [] + ) diff --git a/src/nomad_simulations/schema_packages/particles_state.py b/src/nomad_simulations/schema_packages/particles_state.py index 28afb1e6..4fa5bbe1 100644 --- a/src/nomad_simulations/schema_packages/particles_state.py +++ b/src/nomad_simulations/schema_packages/particles_state.py @@ -5,7 +5,8 @@ import ase.geometry import numpy as np import pint -from deprecated import deprecated + +# from deprecated import deprecated from nomad.datamodel.data import ArchiveSection from nomad.datamodel.metainfo.annotations import ELNAnnotation from nomad.datamodel.metainfo.basesections import Entity @@ -17,6 +18,8 @@ from nomad.metainfo import Context, Section from structlog.stdlib import BoundLogger +from nomad_simulations.schema_packages.atoms_state import State + class Particles: """Particle object. @@ -424,7 +427,7 @@ def _set_positions(self, pos): # ? How generic (usable for any CG model) vs. Martini-specific do we want to be? -class ParticlesState(Entity): +class ParticlesState(State): """ A base section to define individual coarse-grained (CG) particle information. """ @@ -491,3 +494,6 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: # Get particle_type as string, if possible. if not isinstance(self.particle_type, str): self.particle_type = self.resolve_particle_type(logger=logger) + + if self.particle_type is not None: + self.labels = self.particle_type From 12a446f50523bc5481bc125335fc94ad75a5bb26 Mon Sep 17 00:00:00 2001 From: jrudz Date: Thu, 5 Dec 2024 22:40:11 +0100 Subject: [PATCH 2/3] mostly working implementation of particle cell --- .../schema_packages/atoms_state.py | 5 - .../schema_packages/general.py | 29 +- .../schema_packages/model_method.py | 2 +- .../schema_packages/model_system.py | 340 +++++-------- .../schema_packages/particles_state.py | 479 ++---------------- .../schema_packages/utils/utils.py | 1 + tests/test_general.py | 30 +- tests/test_model_method.py | 14 +- 8 files changed, 217 insertions(+), 683 deletions(-) diff --git a/src/nomad_simulations/schema_packages/atoms_state.py b/src/nomad_simulations/schema_packages/atoms_state.py index ed2c6526..a43a2f5a 100644 --- a/src/nomad_simulations/schema_packages/atoms_state.py +++ b/src/nomad_simulations/schema_packages/atoms_state.py @@ -559,7 +559,6 @@ class State(Entity): def __init__(self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs): super().__init__(m_def, m_context, **kwargs) - self.labels = None class AtomsState(State): @@ -651,7 +650,3 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: self.chemical_symbol = self.resolve_chemical_symbol(logger=logger) if self.atomic_number is None: self.atomic_number = self.resolve_atomic_number(logger=logger) - - # Set the labels - if self.chemical_symbol is not None: - self.labels = self.chemical_symbol diff --git a/src/nomad_simulations/schema_packages/general.py b/src/nomad_simulations/schema_packages/general.py index 496cebad..2264e765 100644 --- a/src/nomad_simulations/schema_packages/general.py +++ b/src/nomad_simulations/schema_packages/general.py @@ -225,7 +225,9 @@ def _set_system_branch_depth( ) #! Generalize from checks for atomic systems, error with CG input - def resolve_composition_formula(self, system_parent: ModelSystem) -> None: + def resolve_composition_formula( + self, system_parent: ModelSystem, logger: 'BoundLogger' + ) -> None: """Determine and set the composition formula for `system_parent` and all of its descendants. @@ -246,13 +248,15 @@ def set_composition_formula( to the atom indices stored in system. """ if not subsystems: - atom_indices = ( - system.atom_indices if system.atom_indices is not None else [] + particle_indices = ( + system.particle_indices + if system.particle_indices is not None + else [] ) subsystem_labels = ( - [np.array(labels)[atom_indices]] + [np.array(labels)[particle_indices]] if labels - else ['Unknown' for atom in range(len(atom_indices))] + else ['Unknown' for atom in range(len(particle_indices))] ) else: subsystem_labels = [ @@ -284,18 +288,7 @@ def get_composition_recurs(system: ModelSystem, labels: list[str]) -> None: # ! CG: system_parent.cell[0].particles_state instead of atoms_state! labels = [] if system_parent.cell is not None: - if system_parent.cell[0].name == 'AtomicCell': - labels = ( - [atom.labels for atom in system_parent.cell[0].atoms_state] - if system_parent.cell[0].atoms_state is not None - else [] - ) - elif system_parent.cell[0].name == 'ParticleCell': - labels = ( - [atom.labels for atom in system_parent.cell[0].particles_state] - if system_parent.cell[0].particles_state is not None - else [] - ) + labels = system_parent.cell[0].get('labels', logger=logger) get_composition_recurs(system=system_parent, labels=labels) @@ -322,7 +315,7 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: if is_not_representative(model_system=system_parent, logger=logger): continue - self.resolve_composition_formula(system_parent=system_parent) + self.resolve_composition_formula(system_parent=system_parent, logger=logger) m_package.__init_metainfo__() diff --git a/src/nomad_simulations/schema_packages/model_method.py b/src/nomad_simulations/schema_packages/model_method.py index c7b143ff..7e6cd19d 100644 --- a/src/nomad_simulations/schema_packages/model_method.py +++ b/src/nomad_simulations/schema_packages/model_method.py @@ -544,7 +544,7 @@ def resolve_orbital_references( # If the child is not an "active_atom", the normalization will not run if active_atom.type != 'active_atom': continue - indices = active_atom.atom_indices + indices = active_atom.particle_indices for index in indices: try: active_atoms_state = atoms_state[index] diff --git a/src/nomad_simulations/schema_packages/model_system.py b/src/nomad_simulations/schema_packages/model_system.py index d0cb0042..1c311b66 100644 --- a/src/nomad_simulations/schema_packages/model_system.py +++ b/src/nomad_simulations/schema_packages/model_system.py @@ -22,6 +22,7 @@ from typing import TYPE_CHECKING, Optional import ase +from ase.symbols import symbols2numbers import numpy as np from matid import Classifier, SymmetryAnalyzer # pylint: disable=import-error from matid.classification.classifications import ( @@ -51,7 +52,7 @@ from structlog.stdlib import BoundLogger from nomad_simulations.schema_packages.atoms_state import AtomsState -from nomad_simulations.schema_packages.particles_state import Particles, ParticlesState +from nomad_simulations.schema_packages.particles_state import ParticlesState from nomad_simulations.schema_packages.utils import ( catch_not_implemented, get_sibling_section, @@ -198,6 +199,7 @@ class GeometricSpace(Entity): """, ) + # TODO: Either generalize this or add logic for different cell types def get_geometric_space_for_atomic_cell(self, logger: 'BoundLogger') -> None: """ Get the real space parameters for the atomic cell using ASE. @@ -316,6 +318,8 @@ class Cell(GeometricSpace): """, ) + # ? What does this mean? Number of particles in the cell? + # TODO: improve description n_cell_points = Quantity( type=np.int32, description=""" @@ -328,7 +332,7 @@ class Cell(GeometricSpace): shape=['n_cell_points', 3], unit='meter', description=""" - Positions of all the atoms in Cartesian coordinates. + Positions of all the particles in Cartesian coordinates. """, ) @@ -337,8 +341,8 @@ class Cell(GeometricSpace): shape=['n_cell_points', 3], unit='meter / second', description=""" - Velocities of the atoms. It is the change in cartesian coordinates of the atom position - with time. + Velocities of the particles. It is the change in cartesian coordinates of the + particle position with time. """, ) @@ -347,8 +351,9 @@ class Cell(GeometricSpace): shape=[3, 3], unit='meter', description=""" - Lattice vectors of the simulated cell in Cartesian coordinates. The first index runs - over each lattice vector. The second index runs over the $x, y, z$ Cartesian coordinates. + Lattice vectors of the simulated cell in Cartesian coordinates. The first index + runs over each lattice vector. The second index runs over the $x, y, z$ + Cartesian coordinates. """, ) @@ -371,6 +376,10 @@ class Cell(GeometricSpace): """, ) + def __init__(self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs): + super().__init__(m_def, m_context, **kwargs) + self.logger = None # Initialize logger attribute + @staticmethod def _generate_comparer(obj: 'Cell') -> 'Generator[Any, None, None]': try: @@ -402,6 +411,117 @@ def is_ne_cell(self, other) -> bool: # this does not hold in general, but here we use finite sets return not self.is_equal_cell(other) + def get_state(self): + if isinstance(self, AtomicCell): + return self.atoms_state + elif hasattr(self, 'particles_state'): + return self.particles_state + else: + raise AttributeError( + 'The class does not have atoms_state or particles_state' + ) + + def get(self, key, logger=None): + if key == 'state': + return self.get_state() + elif key == 'labels': + if logger is None: + raise ValueError('Logger is not set') + if isinstance(self, AtomicCell): + return self.get_chemical_symbols(logger) + elif isinstance(self, ParticleCell): + return self.get_particle_types(logger) + else: + return None + else: + raise KeyError(f"Key '{key}' not found in Cell") + + def to_ase_atoms(self, logger: 'BoundLogger') -> 'Optional[ase.Atoms]': + """ + Generates an ASE Atoms object with the most basic information from the parsed `Cell` + section (labels, periodic_boundary_conditions, positions, and lattice_vectors). + + Args: + logger (BoundLogger): The logger to log messages. + + Returns: + (Optional[ase.Atoms]): The ASE Atoms object with the basic information from the `AtomicCell`. + """ + + labels = self.get('labels', logger) + if labels is None: + logger.error('Could not find `Cell.state.labels`.') + labels = ['X'] * len(self.positions) if self.positions is not None else None + # Initialize ase.Atoms object with labels + # ! We need to make sure that the labels from ase.Atoms are not being used downstream! + + try: + symbols2numbers(labels) + except KeyError: + logger.error('Non chemical symbols in `Cell.state.labels`.') + labels = ['X'] * len(labels) + ase_atoms = ase.Atoms(symbols=labels) + + # PBC + if self.periodic_boundary_conditions is None: + logger.info( + 'Could not find `Cell.periodic_boundary_conditions`. They will be set to [False, False, False].' + ) + self.periodic_boundary_conditions = [False, False, False] + ase_atoms.set_pbc(pbc=self.periodic_boundary_conditions) + + # Lattice vectors + if self.lattice_vectors is not None: + ase_atoms.set_cell(cell=self.lattice_vectors.to('angstrom').magnitude) + else: + logger.info('Could not find `AtomicCell.lattice_vectors`.') + + # Positions + if self.positions is not None: + if len(self.positions) != len(self.atoms_state): + logger.error( + 'Length of `Cell.positions` does not coincide with the length of the `Cell.`.' + ) + return None + ase_atoms.set_positions( + newpositions=self.positions.to('angstrom').magnitude + ) + else: + logger.warning('Could not find `AtomicCell.positions`.') + return None + + return ase_atoms + + def from_ase_atoms(self, ase_atoms: ase.Atoms, logger: 'BoundLogger') -> None: + """ + Parses the information from an ASE Atoms object to the `Cell` section. + + Args: + ase_atoms (ase.Atoms): The ASE Atoms object to parse. + logger (BoundLogger): The logger to log messages. + """ + if isinstance(self, AtomicCell): + # `AtomsState[*].chemical_symbol` + for symbol in ase_atoms.get_chemical_symbols(): + atom_state = AtomsState(chemical_symbol=symbol) + self.atoms_state.append(atom_state) + # TODO: implement for `ParticleCell` + + # `periodic_boundary_conditions` + self.periodic_boundary_conditions = ase_atoms.get_pbc() + + # `lattice_vectors` + cell = ase_atoms.get_cell() + self.lattice_vectors = ase.geometry.complete_cell(cell) * ureg('angstrom') + + # `positions` + positions = ase_atoms.get_positions() + if ( + not positions.tolist() + ): # ASE assigns a shape=(0, 3) array if no positions are found + return None + self.positions = positions * ureg('angstrom') + def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: super().normalize(archive, logger) @@ -476,79 +596,6 @@ def get_chemical_symbols(self, logger: 'BoundLogger') -> list[str]: chemical_symbols.append(atom_state.chemical_symbol) return chemical_symbols - def to_ase_atoms(self, logger: 'BoundLogger') -> 'Optional[ase.Atoms]': - """ - Generates an ASE Atoms object with the most basic information from the parsed `AtomicCell` - section (labels, periodic_boundary_conditions, positions, and lattice_vectors). - - Args: - logger (BoundLogger): The logger to log messages. - - Returns: - (Optional[ase.Atoms]): The ASE Atoms object with the basic information from the `AtomicCell`. - """ - # Initialize ase.Atoms object with labels - atoms_labels = self.get_chemical_symbols(logger=logger) - ase_atoms = ase.Atoms(symbols=atoms_labels) - - # PBC - if self.periodic_boundary_conditions is None: - logger.info( - 'Could not find `AtomicCell.periodic_boundary_conditions`. They will be set to [False, False, False].' - ) - self.periodic_boundary_conditions = [False, False, False] - ase_atoms.set_pbc(pbc=self.periodic_boundary_conditions) - - # Lattice vectors - if self.lattice_vectors is not None: - ase_atoms.set_cell(cell=self.lattice_vectors.to('angstrom').magnitude) - else: - logger.info('Could not find `AtomicCell.lattice_vectors`.') - - # Positions - if self.positions is not None: - if len(self.positions) != len(self.atoms_state): - logger.error( - 'Length of `AtomicCell.positions` does not coincide with the length of the `AtomicCell.atoms_state`.' - ) - return None - ase_atoms.set_positions( - newpositions=self.positions.to('angstrom').magnitude - ) - else: - logger.warning('Could not find `AtomicCell.positions`.') - return None - - return ase_atoms - - def from_ase_atoms(self, ase_atoms: ase.Atoms, logger: 'BoundLogger') -> None: - """ - Parses the information from an ASE Atoms object to the `AtomicCell` section. - - Args: - ase_atoms (ase.Atoms): The ASE Atoms object to parse. - logger (BoundLogger): The logger to log messages. - """ - # `AtomsState[*].chemical_symbol` - for symbol in ase_atoms.get_chemical_symbols(): - atom_state = AtomsState(chemical_symbol=symbol) - self.atoms_state.append(atom_state) - - # `periodic_boundary_conditions` - self.periodic_boundary_conditions = ase_atoms.get_pbc() - - # `lattice_vectors` - cell = ase_atoms.get_cell() - self.lattice_vectors = ase.geometry.complete_cell(cell) * ureg('angstrom') - - # `positions` - positions = ase_atoms.get_positions() - if ( - not positions.tolist() - ): # ASE assigns a shape=(0, 3) array if no positions are found - return None - self.positions = positions * ureg('angstrom') - def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: super().normalize(archive, logger) @@ -556,6 +603,7 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: self.name = self.m_def.name if self.name is None else self.name +# TODO Consider changing name to BeadCell or CGBeadCell, only using "particle" for the more abstract reference to atoms or beads class ParticleCell(Cell): """ A base section used to specify the particle cell information of a system. @@ -570,52 +618,11 @@ class ParticleCell(Cell): """, ) - equivalent_particles = Quantity( - type=np.int32, - shape=['n_particle'], - description=""" - List of equivalent particles as defined in `particles`. If no equivalent particles - are found, then the list is simply the index of each element, e.g.: - - [0, 1, 2, 3] all four particles are non-equivalent. - - [0, 0, 0, 3] three equivalent particles and one non-equivalent. - """, - ) - def __init__(self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs): super().__init__(m_def, m_context, **kwargs) # Set the name of the section self.name = self.m_def.name - # def is_equal_cell(self, other) -> bool: - # """ - # Check if the atomic cell is equal to an`other` atomic cell by comparing the `positions` and - # the `AtomsState[*].chemical_symbol`. - # Args: - # other: The other atomic cell to compare with. - # Returns: - # bool: True if the atomic cells are equal, False otherwise. - # """ - # if not isinstance(other, ParticleCell): - # return False - - # # Compare positions using the parent sections's `__eq__` method - # if not super().is_equal_cell(other=other): - # return False - - # # Check that the `chemical_symbol` of the atoms in `cell_1` match with the ones in `cell_2` - # check_positions = self._check_positions( - # positions_1=self.positions, positions_2=other.positions - # ) - # try: - # for particle in check_positions: - # type_1 = self.particles_state[particle[0]].particle_type - # type_2 = other.particles_state[particle[1]].particle_type - # if type_1 != type_2: - # return False - # except Exception: - # return False - # return True - def get_particle_types(self, logger: 'BoundLogger') -> list[str]: """ Get the chemical symbols of the particle in the particle cell. @@ -638,85 +645,6 @@ def get_particle_types(self, logger: 'BoundLogger') -> list[str]: particle_labels.append(particle_state.particle_type) return particle_labels - def to_particles(self, logger: 'BoundLogger') -> Optional[Particles]: - """ - Generates a Particles object with the most basic information from the parsed `ParticleCell` - section (labels, periodic_boundary_conditions, positions, and lattice_vectors). - - Args: - logger (BoundLogger): The logger to log messages. - - Returns: - (Optional[Particles]): The Partilces object with the basic information from the `ParticleCell`. - """ - # Initialize Partilces object with labels - particle_labels = self.get_particle_types(logger=logger) - particles = Particles(symbols=particle_labels) - - # PBC - if self.periodic_boundary_conditions is None: - logger.info( - 'Could not find `ParticleCell.periodic_boundary_conditions`. They will be set to [False, False, False].' - ) - self.periodic_boundary_conditions = [False, False, False] - particles.set_pbc(pbc=self.periodic_boundary_conditions) - - # Lattice vectors - if self.lattice_vectors is not None: - particles.set_cell(cell=self.lattice_vectors.to('angstrom').magnitude) - else: - logger.info('Could not find `ParticleCell.lattice_vectors`.') - - # Positions - if self.positions is not None: - if len(self.positions) != len(self.particles_state): - logger.error( - 'Length of `ParticleCell.positions` does not coincide with the length of the `ParticleCell.particles_state`.' - ) - return None - particles.set_positions( - newpositions=self.positions.to('angstrom').magnitude - ) - else: - logger.warning('Could not find `ParticleCell.positions`.') - return None - - return particles - - def from_particles(self, particles: Particles, logger: 'BoundLogger') -> None: - """ - Parses the information from a Particles object to the `ParticlesCell` section. - - Args: - particles (Particles): The Particles object to parse. - logger (BoundLogger): The logger to log messages. - """ - # `ParticlesState[*].particles_type` - for label in particles.get_particle_types(): - particle_state = ParticlesState(particle_type=label) - self.particles_state.append(particle_state) - - # `periodic_boundary_conditions` - self.periodic_boundary_conditions = particles.get_pbc() - - # `lattice_vectors` - cell = particles.get_cell() - self.lattice_vectors = ase.geometry.complete_cell(cell) * ureg('angstrom') - - # `positions` - positions = particles.get_positions() - if ( - not positions.tolist() - ): # ASE assigns a shape=(0, 3) array if no positions are found - return None - self.positions = positions * ureg('angstrom') - - def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: - super().normalize(archive, logger) - - # Set the name of the section - self.name = self.m_def.name if self.name is None else self.name - class Symmetry(ArchiveSection): """ @@ -1089,7 +1017,7 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: self.m_cache['elemental_composition'] = formula.elemental_composition() -# TODO: generalize! indices instead of atom_indices, state instead of atoms_state... +# TODO work on descriptions to reflect the generalization of ModelSystem class ModelSystem(System): """ Model system used as an input for simulating the material. @@ -1111,7 +1039,7 @@ class ModelSystem(System): formats. This class nest over itself (with the section proxy in `model_system`) to define different - parent-child system trees. The quantities `branch_label`, `branch_depth`, `atom_indices`, + parent-child system trees. The quantities `branch_label`, `branch_depth`, `particle_indices`, and `bond_list` are used to define the parent-child tree. The normalization is ran in the following order: @@ -1232,19 +1160,19 @@ class ModelSystem(System): """, ) - atom_indices = Quantity( + particle_indices = Quantity( type=np.int32, shape=['*'], description=""" - Indices of the atoms in the child with respect to its parent. Example: + Indices of the atoms or, more generally, particles in the child with respect to its parent. Example: - We have SrTiO3, where `AtomicCell.labels = ['Sr', 'Ti', 'O', 'O', 'O']`. If we create a `model_system` child for the `'Ti'` atom only, then in that child - `ModelSystem.model_system.atom_indices = [1]`. If now we want to refer both to - the `'Ti'` and the last `'O'` atoms, `ModelSystem.model_system.atom_indices = [1, 4]`. + `ModelSystem.model_system.particle_indices = [1]`. If now we want to refer both to + the `'Ti'` and the last `'O'` atoms, `ModelSystem.model_system.particle_indices = [1, 4]`. """, ) - # TODO improve description and add an example using the case in atom_indices + # TODO improve description and add an example using the case in particle_indices bond_list = Quantity( type=np.int32, shape=['*', 2], diff --git a/src/nomad_simulations/schema_packages/particles_state.py b/src/nomad_simulations/schema_packages/particles_state.py index 4fa5bbe1..31c0c197 100644 --- a/src/nomad_simulations/schema_packages/particles_state.py +++ b/src/nomad_simulations/schema_packages/particles_state.py @@ -21,411 +21,6 @@ from nomad_simulations.schema_packages.atoms_state import State -class Particles: - """Particle object. - - Adaptation of the ASE Atoms object to coarse-grained particles. For use with - nomad_simulations.model_system.ParticlesCell. - Implemented methods: set_pbc, get_pbc, - set_cell, get_cell, - set_positions, get_positions, - get_particle_types - - Parameters: - - types: str (formula) or list of str. - Default: [‘A’] - typeid: int - (Optional) Id numbers of corresponding particle types. - Default: 0 - type_shapes: str - Store a per-type shape definition for visualization. - A dictionary is stored for each of the NT types, corresponding - to a shape for visualization of that type. - Default: empty - masses: list of float - The mass of each particle. - Default: 1.0 - charges: list of float - Initial atomic charges. - Default: 0.0 - diameter: float - The diameter of each particle. - Default: 1.0 - body: int - The composite body associated with each particle. The value -1 - indicates no body. - Default: -1 - moment_inertia: float - The moment_inertia of each particle (I_xx, I_yy, I_zz). - This inertia tensor is diagonal in the body frame of the particle. - The default value is for point particles. - Default: 0, 0, 0 - positions: float, list of xyz-positions - Particle positions. Needs to be convertible to an - ndarray of shape (N, 3). - Default: 0, 0, 0 - scaled_positions: list of scaled-positions - Like positions, but given in units of the unit cell. - Can not be set at the same time as positions. - Default: 0, 0, 0 - orientation: float - The orientation of each particle. In scalar + vector notation, - this is (r, a_x, a_y, a_z), where the quaternion is q = r + a_xi + a_yj + a_zk. - A unit quaternion has the property: sqrt(r^2 + a_x^2 + a_y^2 + a_z^2) = 1. - Default: 0, 0, 0, 0 - angmom: float - The angular momentum of each particle as a quaternion. - Default: 0, 0, 0, 0 - image: int - The number of times each particle has wrapped around the box (i_x, i_y, i_z). - Default: 0, 0, 0 - cell: 3x3 matrix or length 3 or 6 vector - Unit cell vectors. Can also be given as just three - numbers for orthorhombic cells, or 6 numbers, where - first three are lengths of unit cell vectors, and the - other three are angles between them (in degrees), in following order: - [len(a), len(b), len(c), angle(b,c), angle(a,c), angle(a,b)]. - First vector will lie in x-direction, second in xy-plane, - and the third one in z-positive subspace. - Default value: [0, 0, 0]. - celldisp: Vector - Unit cell displacement vector. To visualize a displaced cell - around the center of mass of a Systems of atoms. Default value - = (0,0,0) - pbc: one or three bool - Periodic boundary conditions flags. Examples: True, - False, 0, 1, (1, 1, 0), (True, False, False). Default - value: False. - - Examples: - - These three are equivalent: - - >>> d = 1.104 # N2 bondlength - >>> a = Atoms('N2', [(0, 0, 0), (0, 0, d)]) - >>> a = Atoms(numbers=[7, 7], positions=[(0, 0, 0), (0, 0, d)]) - >>> a = Atoms([Atom('N', (0, 0, 0)), Atom('N', (0, 0, d))]) - - FCC gold: - - >>> a = 4.05 # Gold lattice constant - >>> b = a / 2 - >>> fcc = Atoms('Au', - ... cell=[(0, b, b), (b, 0, b), (b, b, 0)], - ... pbc=True) - - Hydrogen wire: - - >>> d = 0.9 # H-H distance - >>> h = Atoms('H', positions=[(0, 0, 0)], - ... cell=(d, 0, 0), - ... pbc=(1, 0, 0)) - """ - - def __init__( - self, - types=None, - positions=None, - typeid=None, - type_shapes=None, - moment_inertia=None, - masses=None, - angmom=None, - charges=None, - diameter=None, - body=None, - scaled_positions=None, - orientation=None, - image=None, - cell=None, - pbc=None, - celldisp=None, - ): - self._cellobj = self.set_cell() - self._pbc = np.zeros(3, bool) - - particles = None - - # if hasattr(types, 'get_positions'): - # atoms = types - # types = None - # elif ( - # isinstance(types, (list, tuple)) - # and len(types) > 0 - # and isinstance(types[0], Atom) - # ): - # # Get data from a list or tuple of Atom objects: - # data = [ - # [atom.get_raw(name) for atom in types] - # for name in [ - # 'position', - # 'number', - # 'tag', - # 'momentum', - # 'mass', - # 'magmom', - # 'charge', - # ] - # ] - # atoms = self.__class__(None, *data) - # types = None - - # if atoms is not None: - # # Get data from another Atoms object: - # if scaled_positions is not None: - # raise NotImplementedError - # if types is None and numbers is None: - # numbers = atoms.get_atomic_numbers() - # if positions is None: - # positions = atoms.get_positions() - # if tags is None and atoms.has('tags'): - # tags = atoms.get_tags() - # if momenta is None and atoms.has('momenta'): - # momenta = atoms.get_momenta() - # if magmoms is None and atoms.has('initial_magmoms'): - # magmoms = atoms.get_initial_magnetic_moments() - # if masses is None and atoms.has('masses'): - # masses = atoms.get_masses() - # if charges is None and atoms.has('initial_charges'): - # charges = atoms.get_initial_charges() - # if cell is None: - # cell = atoms.get_cell() - # if celldisp is None: - # celldisp = atoms.get_celldisp() - # if pbc is None: - # pbc = atoms.get_pbc() - - # self.arrays = {} - - # if types is None: - # if numbers is None: - # if positions is not None: - # natoms = len(positions) - # elif scaled_positions is not None: - # natoms = len(scaled_positions) - # else: - # natoms = 0 - # numbers = np.zeros(natoms, int) - # self.new_array('numbers', numbers, int) - # else: - # if numbers is not None: - # raise TypeError('Use only one of "types" and "numbers".') - # else: - # self.new_array('numbers', types2numbers(types), int) - - # if self.numbers.ndim != 1: - # raise ValueError('"numbers" must be 1-dimensional.') - - if cell is None: - cell = np.zeros((3, 3)) - self.set_cell(cell) - - # if celldisp is None: - # celldisp = np.zeros(shape=(3, 1)) - # self.set_celldisp(celldisp) - - # if positions is None: - # if scaled_positions is None: - # positions = np.zeros((len(self.arrays['numbers']), 3)) - # else: - # assert self.cell.rank == 3 - # positions = np.dot(scaled_positions, self.cell) - # else: - # if scaled_positions is not None: - # raise TypeError('Use only one of "types" and "numbers".') - # self.new_array('positions', positions, float, (3,)) - # self.set_tags(default(tags, 0)) - # self.set_masses(default(masses, None)) - # self.set_initial_magnetic_moments(default(magmoms, 0.0)) - # self.set_initial_charges(default(charges, 0.0)) - # if pbc is None: - # pbc = False - # self.set_pbc(pbc) - # self.set_momenta(default(momenta, (0.0, 0.0, 0.0)), apply_constraint=False) - - # if velocities is not None: - # if momenta is None: - # self.set_velocities(velocities) - # else: - # raise TypeError('Use only one of "momenta" and "velocities".') - - # if info is None: - # self.info = {} - # else: - # self.info = dict(info) - - # self.calc = calculator - - # def set_cell(self, cell): - # if cell is None: - # cell = np.zeros((3, 3)) - - # @property - # def symbols(self): - # """Get chemical symbols as a :class:`ase.symbols.Symbols` object. - - # The object works like ``atoms.numbers`` except its values - # are strings. It supports in-place editing.""" - # return Symbols(self.numbers) - - # @symbols.setter - # def symbols(self, obj): - # new_symbols = Symbols.fromsymbols(obj) - # self.numbers[:] = new_symbols.numbers - - def get_particle_types(self): - """Get list of particle type strings. - - Labels describing type of coarse-grained particles.""" - return list(self.types) - - def set_cell(self, cell, scale_atoms=False, apply_constraint=True): - """Set unit cell vectors. - - Parameters: - - cell: 3x3 matrix - Unit cell. A 3x3 matrix (the three unit cell vectors). - First vector will lie in x-direction, second in - xy-plane, and the third one in z-positive subspace. - scale_atoms: bool - Fix atomic positions or move atoms with the unit cell? - Default behavior is to *not* move the atoms (scale_atoms=False). - apply_constraint: bool - Whether to apply constraints to the given cell. - """ - - # Override pbcs if and only if given a Cell object: - cell = ase.Cell.new(cell) - - if apply_constraint and hasattr(self, '_constraints'): - for constraint in self.constraints: - if hasattr(constraint, 'adjust_cell'): - constraint.adjust_cell(self, cell) - - if scale_atoms: - M = np.linalg.solve(self.cell.complete(), cell.complete()) - self.positions[:] = np.dot(self.positions, M) - - self.cell[:] = cell - - def get_cell(self, complete=False): - """Get the three unit cell vectors as a `class`:ase.cell.Cell` object. - - The Cell object resembles a 3x3 ndarray, and cell[i, j] - is the jth Cartesian coordinate of the ith cell vector.""" - if complete: - cell = self.cell.complete() - else: - cell = self.cell.copy() - - return cell - - @property - def pbc(self): - """Reference to pbc-flags for in-place manipulations.""" - return self._pbc - - @pbc.setter - def pbc(self, pbc): - self._pbc[:] = pbc - - def set_pbc(self, pbc): - """Set periodic boundary condition flags.""" - self.pbc = pbc - - def get_pbc(self): - """Get periodic boundary condition flags.""" - return self.pbc.copy() - - def set_positions(self, newpositions, apply_constraint=True): - """Set positions, honoring any constraints. To ignore constraints, - use *apply_constraint=False*.""" - if self.constraints and apply_constraint: - newpositions = np.array(newpositions, float) - for constraint in self.constraints: - constraint.adjust_positions(self, newpositions) - - self.set_array('positions', newpositions, shape=(3,)) - - def get_positions(self, wrap=False, **wrap_kw): - """Get array of positions. - - Parameters: - - wrap: bool - wrap atoms back to the cell before returning positions - wrap_kw: (keyword=value) pairs - optional keywords `pbc`, `center`, `pretty_translation`, `eps`, - see :func:`ase.geometry.wrap_positions` - """ - if wrap: - if 'pbc' not in wrap_kw: - wrap_kw['pbc'] = self.pbc - return ase.geometry.wrap_positions(self.positions, self.cell, **wrap_kw) - else: - return self.arrays['positions'].copy() - - def get_scaled_positions(self, wrap=True): - """Get positions relative to unit cell. - - If wrap is True, atoms outside the unit cell will be wrapped into - the cell in those directions with periodic boundary conditions - so that the scaled coordinates are between zero and one. - - If any cell vectors are zero, the corresponding coordinates - are evaluated as if the cell were completed using - ``cell.complete()``. This means coordinates will be Cartesian - as long as the non-zero cell vectors span a Cartesian axis or - plane.""" - - fractional = self.cell.scaled_positions(self.positions) - - if wrap: - for i, periodic in enumerate(self.pbc): - if periodic: - # Yes, we need to do it twice. - # See the scaled_positions.py test. - fractional[:, i] %= 1.0 - fractional[:, i] %= 1.0 - - return fractional - - def set_scaled_positions(self, scaled): - """Set positions relative to unit cell.""" - self.positions[:] = self.cell.cartesian_positions(scaled) - - def wrap(self, **wrap_kw): - """Wrap positions to unit cell. - - Parameters: - - wrap_kw: (keyword=value) pairs - optional keywords `pbc`, `center`, `pretty_translation`, `eps`, - see :func:`ase.geometry.wrap_positions` - """ - - if 'pbc' not in wrap_kw: - wrap_kw['pbc'] = self.pbc - - self.positions[:] = self.get_positions(wrap=True, **wrap_kw) - - def _get_positions(self): - """Return reference to positions-array for in-place manipulations.""" - return self.arrays['positions'] - - def _set_positions(self, pos): - """Set positions directly, bypassing constraints.""" - self.arrays['positions'][:] = pos - - positions = property( - _get_positions, - _set_positions, - doc='Attribute for direct ' + 'manipulation of the positions.', - ) - - # ? How generic (usable for any CG model) vs. Martini-specific do we want to be? class ParticlesState(State): """ @@ -436,40 +31,65 @@ class ParticlesState(State): particle_type = Quantity( type=str, description=""" - Symbol(s) describing the CG particle type. Currently, entrie particle label is + Symbol(s) describing the CG particle type. Currently, entire particle label is used for type definition. """, ) - # ? Do we want to reflect the Martini size nomenclature and include bead volume/bead mass? - # particle_size = Quantity( - # type=np.float64, - # description=""" - # Particle size, determining the number of non-hydrogen atoms represented by the - # particle. Currently, possible values are 0.47 nm (regular, default), - # 0.43/0.41 nm (small), and 0.34 nm (tiny). - # """, - # ) + mass = Quantity( + type=np.float64, + unit='kg', + description=""" + Total mass of the particle. + """, + ) - # particle_mass = Quantity( - # type=np.float64, - # description=""" - # Particle size, determining the number of non-hydrogen atoms represented by the - # particle. Currently, possible values are 72 amu (regular, default), 54/45 amu - # (small), and 36 amu (tiny). - # """, - # ) + charge = Quantity( + type=np.float64, + unit='coulomb', + description=""" + Total charge of the particle. + """, + ) charge = Quantity( - type=np.int32, - default=0, + type=np.float64, + unit='coulomb', description=""" - Charge of the particle. Neutral = 0. Can be any positive integer (+1, +2...) - for cations or any negative integer (-1, -2...) for anions. + Total charge of the particle. """, - a_eln=ELNAnnotation(component='NumberEditQuantity'), ) + # Other possible quantities + # diameter: float + # The diameter of each particle. + # Default: 1.0 + # body: int + # The composite body associated with each particle. The value -1 + # indicates no body. + # Default: -1 + # moment_inertia: float + # The moment_inertia of each particle (I_xx, I_yy, I_zz). + # This inertia tensor is diagonal in the body frame of the particle. + # The default value is for point particles. + # Default: 0, 0, 0 + # scaled_positions: list of scaled-positions #! for cell if relevant + # Like positions, but given in units of the unit cell. + # Can not be set at the same time as positions. + # Default: 0, 0, 0 + # orientation: float + # The orientation of each particle. In scalar + vector notation, + # this is (r, a_x, a_y, a_z), where the quaternion is q = r + a_xi + a_yj + a_zk. + # A unit quaternion has the property: sqrt(r^2 + a_x^2 + a_y^2 + a_z^2) = 1. + # Default: 0, 0, 0, 0 + # angmom: float #? for cell or here? + # The angular momentum of each particle as a quaternion. + # Default: 0, 0, 0, 0 + # image: int #! advance PBC stuff would go in cell I guess + # The number of times each particle has wrapped around the box (i_x, i_y, i_z). + # Default: 0, 0, 0 + + # ? What is the purpose exactly of this function? Example? def resolve_particle_type(self, logger: 'BoundLogger') -> Optional[str]: """ Checks if any value is passed as particle label. Converts to string to be used as @@ -494,6 +114,3 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None: # Get particle_type as string, if possible. if not isinstance(self.particle_type, str): self.particle_type = self.resolve_particle_type(logger=logger) - - if self.particle_type is not None: - self.labels = self.particle_type diff --git a/src/nomad_simulations/schema_packages/utils/utils.py b/src/nomad_simulations/schema_packages/utils/utils.py index eff18376..6483c428 100644 --- a/src/nomad_simulations/schema_packages/utils/utils.py +++ b/src/nomad_simulations/schema_packages/utils/utils.py @@ -48,6 +48,7 @@ def get_sibling_section( if not sibling_section_name: logger.warning('The sibling_section_name is empty.') return None + sibling_section = section.m_xpath(f'm_parent.{sibling_section_name}', dict=False) # If the sibling_section is a list, return the element `index_sibling` of that list if isinstance(sibling_section, list): diff --git a/tests/test_general.py b/tests/test_general.py index a693f33c..621db4c9 100644 --- a/tests/test_general.py +++ b/tests/test_general.py @@ -87,7 +87,7 @@ def get_flat_depths( assert value == result @pytest.mark.parametrize( - 'is_representative, has_atom_indices, mol_label_list, n_mol_list, atom_labels_list, composition_formula_list, custom_formulas', + 'is_representative, has_particle_indices, mol_label_list, n_mol_list, atom_labels_list, composition_formula_list, custom_formulas', [ ( True, @@ -168,7 +168,7 @@ def get_flat_depths( def test_system_hierarchy_for_molecules( self, is_representative: bool, - has_atom_indices: bool, + has_particle_indices: bool, mol_label_list: list[str], n_mol_list: list[int], atom_labels_list: list[str], @@ -181,8 +181,8 @@ def test_system_hierarchy_for_molecules( Args: is_representative (bool): Specifies if branch_depth = 0 is representative or not. If not representative, the composition formulas should not be generated. - has_atom_indices (bool): Specifies if the atom_indices should be populated during parsing. - Without atom_indices, the composition formulas for the deepest level of the hierarchy + has_particle_indices (bool): Specifies if the particle_indices should be populated during parsing. + Without particle_indices, the composition formulas for the deepest level of the hierarchy should not be populated. mol_label_list (list[str]): Molecule types for generating the hierarchy. n_mol_list (list[int]): Number of molecules for each molecule type. Should be same @@ -212,15 +212,15 @@ def test_system_hierarchy_for_molecules( ctr_comp = 1 atomic_cell = AtomicCell() model_system.cell.append(atomic_cell) - if has_atom_indices: - model_system.atom_indices = [] + if has_particle_indices: + model_system.particle_indices = [] for mol_label, n_mol, atom_labels in zip( mol_label_list, n_mol_list, atom_labels_list ): # Create a branch in the hierarchy for this molecule type model_system_mol_group = ModelSystem() - if has_atom_indices: - model_system_mol_group.atom_indices = [] + if has_particle_indices: + model_system_mol_group.particle_indices = [] model_system_mol_group.branch_label = ( f'group_{mol_label}' if mol_label is not None else None ) @@ -241,14 +241,14 @@ def test_system_hierarchy_for_molecules( AtomsState(chemical_symbol=atom_label) ) n_atoms = len(atomic_cell.atoms_state) - atom_indices = np.arange(n_atoms - len(atom_labels), n_atoms) - if has_atom_indices: - model_system_mol.atom_indices = atom_indices - model_system_mol_group.atom_indices = np.append( - model_system_mol_group.atom_indices, atom_indices + particle_indices = np.arange(n_atoms - len(atom_labels), n_atoms) + if has_particle_indices: + model_system_mol.particle_indices = particle_indices + model_system_mol_group.particle_indices = np.append( + model_system_mol_group.particle_indices, particle_indices ) - model_system.atom_indices = np.append( - model_system.atom_indices, atom_indices + model_system.particle_indices = np.append( + model_system.particle_indices, particle_indices ) simulation.normalize(EntryArchive(), logger) diff --git a/tests/test_model_method.py b/tests/test_model_method.py index 6f6e6393..d352b858 100644 --- a/tests/test_model_method.py +++ b/tests/test_model_method.py @@ -79,7 +79,7 @@ def test_resolve_type(self, tb_section: TB, result: Optional[str]): is_representative=True, cell=[AtomicCell(atoms_state=[AtomsState()])], model_system=[ - ModelSystem(type='active_atom', atom_indices=[2]) + ModelSystem(type='active_atom', particle_indices=[2]) ], ) ], @@ -93,7 +93,7 @@ def test_resolve_type(self, tb_section: TB, result: Optional[str]): is_representative=True, cell=[AtomicCell(atoms_state=[AtomsState(orbitals_state=[])])], model_system=[ - ModelSystem(type='active_atom', atom_indices=[0]) + ModelSystem(type='active_atom', particle_indices=[0]) ], ) ], @@ -117,7 +117,7 @@ def test_resolve_type(self, tb_section: TB, result: Optional[str]): ) ], model_system=[ - ModelSystem(type='active_atom', atom_indices=[0]) + ModelSystem(type='active_atom', particle_indices=[0]) ], ) ], @@ -207,7 +207,7 @@ def test_resolve_orbital_references( is_representative=True, cell=[AtomicCell(atoms_state=[AtomsState()])], model_system=[ - ModelSystem(type='active_atom', atom_indices=[2]) + ModelSystem(type='active_atom', particle_indices=[2]) ], ) ], @@ -222,7 +222,7 @@ def test_resolve_orbital_references( is_representative=True, cell=[AtomicCell(atoms_state=[AtomsState(orbitals_state=[])])], model_system=[ - ModelSystem(type='active_atom', atom_indices=[0]) + ModelSystem(type='active_atom', particle_indices=[0]) ], ) ], @@ -247,7 +247,7 @@ def test_resolve_orbital_references( ) ], model_system=[ - ModelSystem(type='active_atom', atom_indices=[0]) + ModelSystem(type='active_atom', particle_indices=[0]) ], ) ], @@ -272,7 +272,7 @@ def test_resolve_orbital_references( ) ], model_system=[ - ModelSystem(type='active_atom', atom_indices=[0]) + ModelSystem(type='active_atom', particle_indices=[0]) ], ) ], From 163f846f500ddc3997c0a71f0733721e46d3144c Mon Sep 17 00:00:00 2001 From: jrudz Date: Tue, 10 Dec 2024 09:32:38 +0100 Subject: [PATCH 3/3] errors to warnings --- .../schema_packages/model_system.py | 24 +++++++++++++------ 1 file changed, 17 insertions(+), 7 deletions(-) diff --git a/src/nomad_simulations/schema_packages/model_system.py b/src/nomad_simulations/schema_packages/model_system.py index 1c311b66..eecc1afb 100644 --- a/src/nomad_simulations/schema_packages/model_system.py +++ b/src/nomad_simulations/schema_packages/model_system.py @@ -450,16 +450,26 @@ def to_ase_atoms(self, logger: 'BoundLogger') -> 'Optional[ase.Atoms]': labels = self.get('labels', logger) if labels is None: - logger.error('Could not find `Cell.state.labels`.') + # ! Check the scope of this message + logger.error( + 'Could not find `Cell.state.labels`.' + 'Using ase functionalities with `X` as labels.' + 'This is normal for non-atomic particles,' + 'but no particle labels will be stored.' + ) labels = ['X'] * len(self.positions) if self.positions is not None else None # Initialize ase.Atoms object with labels # ! We need to make sure that the labels from ase.Atoms are not being used downstream! - - try: - symbols2numbers(labels) - except KeyError: - logger.error('Non chemical symbols in `Cell.state.labels`.') - labels = ['X'] * len(labels) + else: + try: + symbols2numbers(labels) + except KeyError: + logger.warning( + 'Non chemical symbols in `Cell.state.labels`.' + 'Using ase functionalities with `X` as labels.' + 'This is normal for non-atomic particles.' + ) + labels = ['X'] * len(labels) ase_atoms = ase.Atoms(symbols=labels) # PBC