From 6ca7d50f0346e8133b8d7ab717df38d6b3994dec Mon Sep 17 00:00:00 2001 From: chrisjonesbsu Date: Wed, 8 May 2024 14:40:31 -0600 Subject: [PATCH 1/3] add ruff changes --- .github/workflows/pytest.yml | 2 +- .pre-commit-config.yaml | 29 +-- docs/source/conf.py | 4 +- flowermd/base/forcefield.py | 8 +- flowermd/base/molecule.py | 29 +-- flowermd/base/simulation.py | 47 +--- flowermd/base/system.py | 82 ++----- flowermd/internal/ff_utils.py | 4 +- flowermd/internal/utils.py | 4 +- flowermd/library/forcefields.py | 29 +-- flowermd/library/polymers.py | 12 +- flowermd/library/simulations/tensile.py | 8 +- .../surface_wetting/surface_wetting.py | 25 +-- flowermd/modules/surface_wetting/utils.py | 4 +- flowermd/modules/welding/utils.py | 8 +- flowermd/modules/welding/welding.py | 24 +- flowermd/tests/base/test_molecule.py | 40 +--- flowermd/tests/base/test_simulation.py | 20 +- flowermd/tests/base/test_system.py | 210 +++++------------- flowermd/tests/base_test.py | 18 +- flowermd/tests/library/test_forcefields.py | 16 +- flowermd/tests/library/test_tensile.py | 4 +- .../tests/modules/test_surface_wetting.py | 19 +- flowermd/tests/modules/test_weld.py | 16 +- flowermd/tests/utils/test_actions.py | 12 +- flowermd/tests/utils/test_rigid_body.py | 3 +- flowermd/tests/utils/test_utils.py | 28 +-- flowermd/utils/rigid_body.py | 16 +- flowermd/utils/utils.py | 4 +- setup.py | 8 +- 30 files changed, 193 insertions(+), 540 deletions(-) diff --git a/.github/workflows/pytest.yml b/.github/workflows/pytest.yml index 3c08cfb4..284423a0 100644 --- a/.github/workflows/pytest.yml +++ b/.github/workflows/pytest.yml @@ -2,7 +2,6 @@ name: pytest on: push: - # Action will run when any changes to these paths are pushed or pr'ed to master branches: [ main ] paths: - flowermd/** @@ -26,6 +25,7 @@ on: jobs: pytest: + if: github.event.pull_request.draft == false strategy: fail-fast: false matrix: diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 63a4ad32..c6c5daf3 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -1,15 +1,19 @@ ci: autofix_commit_msg: | [pre-commit.ci] auto fixes from pre-commit.com hooks - for more information, see https://pre-commit.ci autofix_prs: true autoupdate_commit_msg: '[pre-commit.ci] pre-commit autoupdate' autoupdate_schedule: weekly skip: [ ] submodules: false - repos: + - repo: https://github.com/astral-sh/ruff-pre-commit + rev: v0.4.2 # Ruff version + hooks: + - id: ruff + args: [--line-length=80, --fix, --extend-ignore=E203] + - id: ruff-format - repo: https://github.com/pre-commit/pre-commit-hooks rev: v4.6.0 hooks: @@ -17,11 +21,6 @@ repos: - id: end-of-file-fixer - id: trailing-whitespace exclude: 'flowermd/tests/assets/.* | flowermd/assets/.*' - - repo: https://github.com/psf/black - rev: 24.4.2 - hooks: - - id: black - args: [ --line-length=80 ] - repo: https://github.com/pycqa/isort rev: 5.13.2 hooks: @@ -30,19 +29,3 @@ repos: args: [ --profile=black, --line-length=80 ] exclude: 'flowermd/tests/assets/.* ' - - - repo: https://github.com/pycqa/flake8 - rev: 7.0.0 - hooks: - - id: flake8 - args: - - --max-line-length=80 - - --extend-ignore=E203 - exclude: '__init__.py' - - - repo: https://github.com/pycqa/pydocstyle - rev: '6.3.0' - hooks: - - id: pydocstyle - exclude: ^(flowermd/tests/|flowermd/internal/|flowermd/utils|setup.py|flowermd/__version__.py|docs/) - args: [ --convention=numpy ] diff --git a/docs/source/conf.py b/docs/source/conf.py index 18e21b53..a6ccd74a 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -10,9 +10,7 @@ import sys project = "flowerMD" -copyright = ( - "2023, Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" -) +copyright = "2023, Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" author = "Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" sys.path.insert(0, os.path.abspath("../..")) diff --git a/flowermd/base/forcefield.py b/flowermd/base/forcefield.py index 532669d4..62dff440 100644 --- a/flowermd/base/forcefield.py +++ b/flowermd/base/forcefield.py @@ -11,9 +11,7 @@ def __init__(self, forcefield_files=None, name=None): super(BaseXMLForcefield, self).__init__( forcefield_files=forcefield_files, name=name ) - self.gmso_ff = ( - ffutils.FoyerFFs().load(forcefield_files or name).to_gmso_ff() - ) + self.gmso_ff = ffutils.FoyerFFs().load(forcefield_files or name).to_gmso_ff() class BaseHOOMDForcefield: @@ -22,8 +20,6 @@ class BaseHOOMDForcefield: def __init__(self, hoomd_forces): self.hoomd_forces = hoomd_forces if hoomd_forces is None: - raise NotImplementedError( - "`hoomd_forces` must be defined in the subclass." - ) + raise NotImplementedError("`hoomd_forces` must be defined in the subclass.") if not isinstance(hoomd_forces, list): raise TypeError("`hoomd_forces` must be a list.") diff --git a/flowermd/base/molecule.py b/flowermd/base/molecule.py index 840d32b9..a37d2a3c 100644 --- a/flowermd/base/molecule.py +++ b/flowermd/base/molecule.py @@ -204,8 +204,7 @@ def _load(self): return mb.load(self.smiles, smiles=True) else: raise MoleculeLoadError( - msg=f"Unable to load the molecule from smiles " - f"{self.smiles}." + msg=f"Unable to load the molecule from smiles " f"{self.smiles}." ) def _align_backbones_z_axis(self, heavy_atoms_only=False): @@ -214,11 +213,7 @@ def _align_backbones_z_axis(self, heavy_atoms_only=False): if heavy_atoms_only: try: positions = np.array( - [ - p.xyz[0] - for p in mol.particles() - if p.element.symbol != "H" - ] + [p.xyz[0] for p in mol.particles() if p.element.symbol != "H"] ) except AttributeError: positions = mol.xyz @@ -277,9 +272,7 @@ def _identify_particle_information(self, gmso_molecule): ): self.hydrogen_types.append(p_name) self.particle_typeid.append(self.particle_types.index(p_name)) - self.particle_charge.append( - site.charge.to_value() if site.charge else 0 - ) + self.particle_charge.append(site.charge.to_value() if site.charge else 0) def _identify_pairs(self, particle_types): """Identify all unique particle pairs from the particle types. @@ -290,9 +283,7 @@ def _identify_pairs(self, particle_types): List of all particle types. """ - self.pairs = set( - itertools.combinations_with_replacement(particle_types, 2) - ) + self.pairs = set(itertools.combinations_with_replacement(particle_types, 2)) def _identify_bond_types(self, gmso_molecule): """Identify all unique bond types from the GMSO topology. @@ -314,7 +305,7 @@ def _identify_bond_types(self, gmso_molecule): or bond.connection_members[1].name ) bond_connections = [p1_name, p2_name] - if not tuple(bond_connections[::-1]) in self.bond_types: + if tuple(bond_connections[::-1]) not in self.bond_types: self.bond_types.add(tuple(bond_connections)) def _identify_angle_types(self, gmso_molecule): @@ -341,7 +332,7 @@ def _identify_angle_types(self, gmso_molecule): or angle.connection_members[2].name ) angle_connections = [p1_name, p2_name, p3_name] - if not tuple(angle_connections[::-1]) in self.angle_types: + if tuple(angle_connections[::-1]) not in self.angle_types: self.angle_types.add(tuple(angle_connections)) def _identify_dihedral_types(self, gmso_molecule): @@ -372,7 +363,7 @@ def _identify_dihedral_types(self, gmso_molecule): or dihedral.connection_members[3].name ) dihedral_connections = [p1_name, p2_name, p3_name, p4_name] - if not tuple(dihedral_connections[::-1]) in self.dihedral_types: + if tuple(dihedral_connections[::-1]) not in self.dihedral_types: self.dihedral_types.add(tuple(dihedral_connections)) def _identify_improper_types(self, gmso_molecule): @@ -403,7 +394,7 @@ def _identify_improper_types(self, gmso_molecule): or improper.connection_members[3].name ) improper_connections = [p1_name, p2_name, p3_name, p4_name] - if not tuple(improper_connections[::-1]) in self.improper_types: + if tuple(improper_connections[::-1]) not in self.improper_types: self.improper_types.add(tuple(improper_connections)) def _identify_topology_information(self, gmso_molecule): @@ -435,9 +426,7 @@ def _validate_force_field(self): # Update topology information from typed gmso after applying ff. self._identify_topology_information(self.gmso_molecule) elif isinstance(self.force_field, BaseHOOMDForcefield): - _validate_hoomd_ff( - self.force_field.hoomd_forces, self.topology_information - ) + _validate_hoomd_ff(self.force_field.hoomd_forces, self.topology_information) elif isinstance(self.force_field, List): _validate_hoomd_ff(self.force_field, self.topology_information) else: diff --git a/flowermd/base/simulation.py b/flowermd/base/simulation.py index 04ba2d44..2afcbc4d 100644 --- a/flowermd/base/simulation.py +++ b/flowermd/base/simulation.py @@ -71,8 +71,7 @@ def __init__( ): if not isinstance(forcefield, Iterable) or isinstance(forcefield, str): raise ValueError( - "forcefield must be a sequence of " - "hoomd.md.force.Force objects." + "forcefield must be a sequence of " "hoomd.md.force.Force objects." ) else: for obj in forcefield: @@ -87,9 +86,7 @@ def __init__( self.gsd_write_freq = int(gsd_write_freq) self.maximum_write_buffer_size = gsd_max_buffer_size self.log_write_freq = int(log_write_freq) - self._std_out_freq = int( - (self.gsd_write_freq + self.log_write_freq) / 2 - ) + self._std_out_freq = int((self.gsd_write_freq + self.log_write_freq) / 2) self.gsd_file_name = gsd_file_name self.log_file_name = log_file_name self.log_quantities = [ @@ -436,9 +433,7 @@ def thermostat(self, thermostat): thermostat : flowermd.utils.HOOMDThermostats, required The type of thermostat to use. """ - if not issubclass( - self._thermostat, hoomd.md.methods.thermostats.Thermostat - ): + if not issubclass(self._thermostat, hoomd.md.methods.thermostats.Thermostat): raise ValueError( f"Invalid thermostat. Please choose from: {HOOMDThermostats}" ) @@ -527,9 +522,7 @@ def _initialize_thermostat(self, thermostat_kwargs): required_thermostat_kwargs = {} for k in inspect.signature(self.thermostat).parameters: if k not in thermostat_kwargs.keys(): - raise ValueError( - f"Missing required parameter {k} for thermostat." - ) + raise ValueError(f"Missing required parameter {k} for thermostat.") required_thermostat_kwargs[k] = thermostat_kwargs[k] return self.thermostat(**required_thermostat_kwargs) @@ -551,9 +544,7 @@ def set_integrator_method(self, integrator_method, method_kwargs): if not self.integrator: # Integrator and method not yet created self.integrator = hoomd.md.Integrator( dt=self.dt, - integrate_rotational_dof=( - True if self._rigid_constraint else False - ), + integrate_rotational_dof=(True if self._rigid_constraint else False), ) if self._rigid_constraint: self.integrator.rigid = self._rigid_constraint @@ -715,9 +706,7 @@ def run_update_volume( self.set_integrator_method( integrator_method=hoomd.md.methods.ConstantVolume, method_kwargs={ - "thermostat": self._initialize_thermostat( - {"kT": kT, "tau": tau_kt} - ), + "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), "filter": self.integrate_group, }, ) @@ -849,9 +838,7 @@ def run_NPT( "rescale_all": rescale_all, "gamma": gamma, "filter": self.integrate_group, - "thermostat": self._initialize_thermostat( - {"kT": kT, "tau": tau_kt} - ), + "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), }, ) if thermalize_particles: @@ -894,9 +881,7 @@ def run_NVT( self.set_integrator_method( integrator_method=hoomd.md.methods.ConstantVolume, method_kwargs={ - "thermostat": self._initialize_thermostat( - {"kT": kT, "tau": tau_kt} - ), + "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), "filter": self.integrate_group, }, ) @@ -1105,17 +1090,13 @@ def _thermalize_system(self, kT): filter=self.integrate_group, kT=kT.range[0] ) else: - self.state.thermalize_particle_momenta( - filter=self.integrate_group, kT=kT - ) + self.state.thermalize_particle_momenta(filter=self.integrate_group, kT=kT) def _lj_force(self): """Return the Lennard-Jones pair force.""" if not self.integrator: lj_force = [ - f - for f in self._forcefield - if isinstance(f, hoomd.md.pair.pair.LJ) + f for f in self._forcefield if isinstance(f, hoomd.md.pair.pair.LJ) ][0] else: lj_force = [ @@ -1141,9 +1122,7 @@ def _create_state(self, initial_state): print("Initializing simulation state from a GSD file.") self.create_state_from_gsd(initial_state) elif isinstance(initial_state, hoomd.snapshot.Snapshot): - print( - "Initializing simulation state from a hoomd.snapshot.Snapshot" - ) + print("Initializing simulation state from a hoomd.snapshot.Snapshot") self.create_state_from_snapshot(initial_state) elif isinstance(initial_state, gsd.hoomd.Frame): print("Initializing simulation state from a gsd.hoomd.Frame.") @@ -1151,9 +1130,7 @@ def _create_state(self, initial_state): def _add_hoomd_writers(self): """Create gsd and log writers.""" - gsd_logger = hoomd.logging.Logger( - categories=["scalar", "string", "sequence"] - ) + gsd_logger = hoomd.logging.Logger(categories=["scalar", "string", "sequence"]) logger = hoomd.logging.Logger(categories=["scalar", "string"]) gsd_logger.add(self, quantities=["timestep", "tps"]) logger.add(self, quantities=["timestep", "tps"]) diff --git a/flowermd/base/system.py b/flowermd/base/system.py index 6416be53..a9609160 100644 --- a/flowermd/base/system.py +++ b/flowermd/base/system.py @@ -89,16 +89,12 @@ def __init__( if isinstance(mol_item, Molecule): # keep track of molecule types indices to assign to sites # before applying forcefield - self._mol_type_idx.extend( - [self.n_mol_types] * mol_item.n_particles - ) + self._mol_type_idx.extend([self.n_mol_types] * mol_item.n_particles) self.all_molecules.extend(mol_item.molecules) # if ff is provided in the Molecule class use that as the ff if mol_item.force_field: if isinstance(mol_item.force_field, BaseHOOMDForcefield): - self._hoomd_forcefield.extend( - mol_item.force_field.hoomd_forces - ) + self._hoomd_forcefield.extend(mol_item.force_field.hoomd_forces) elif isinstance(mol_item.force_field, BaseXMLForcefield): self._gmso_forcefields_dict[str(self.n_mol_types)] = ( mol_item.force_field.gmso_ff @@ -155,9 +151,7 @@ def mass(self): @property def net_charge(self): """Net charge of the system.""" - return sum( - site.charge if site.charge else 0 for site in self.gmso_system.sites - ) + return sum(site.charge if site.charge else 0 for site in self.gmso_system.sites) @property def box(self): @@ -309,13 +303,8 @@ def hoomd_snapshot(self): @property def hoomd_forcefield(self): """List of HOOMD forces.""" - if ( - self._ff_refs != self.reference_values - and self._gmso_forcefields_dict - ): - self._hoomd_forcefield = self._create_hoomd_forcefield( - **self._ff_kwargs - ) + if self._ff_refs != self.reference_values and self._gmso_forcefields_dict: + self._hoomd_forcefield = self._create_hoomd_forcefield(**self._ff_kwargs) return self._hoomd_forcefield def remove_hydrogens(self): @@ -327,9 +316,7 @@ def remove_hydrogens(self): """ # Try by element first: hydrogens = [ - site - for site in self.gmso_system.sites - if site.element.atomic_number == 1 + site for site in self.gmso_system.sites if site.element.atomic_number == 1 ] # If none found by element; try by mass if len(hydrogens) == 0: @@ -339,9 +326,7 @@ def remove_hydrogens(self): if site.mass.to("amu").value == 1.008 ] if len(hydrogens) == 0: - warnings.warn( - "Hydrogen atoms could not be found by element or mass" - ) + warnings.warn("Hydrogen atoms could not be found by element or mass") for h in hydrogens: # Find bond and other site in bond, add mass and charge for bond in self.gmso_system.iter_connections_by_site( @@ -365,18 +350,15 @@ def _scale_charges(self): """ charges = np.array( - [ - site.charge if site.charge else 0 - for site in self.gmso_system.sites - ] + [site.charge if site.charge else 0 for site in self.gmso_system.sites] ) net_charge = sum(charges) abs_charge = sum(abs(charges)) if abs_charge != 0: for site in self.gmso_system.sites: - site.charge -= abs( - site.charge if site.charge else 0 * u.Unit("C") - ) * (net_charge / abs_charge) + site.charge -= abs(site.charge if site.charge else 0 * u.Unit("C")) * ( + net_charge / abs_charge + ) def pickle_forcefield(self, file_path="forcefield.pickle"): """Pickle the list of HOOMD forces. @@ -419,9 +401,7 @@ def _create_hoomd_forcefield(self, r_cut, nlist_buffer, pppm_kwargs): nlist_buffer=nlist_buffer, pppm_kwargs=pppm_kwargs, auto_scale=False, - base_units=( - self._reference_values if self._reference_values else None - ), + base_units=(self._reference_values if self._reference_values else None), ) for force in ff: force_list.extend(ff[force]) @@ -433,9 +413,7 @@ def _create_hoomd_snapshot(self): snap, refs = to_gsd_snapshot( top=self.gmso_system, auto_scale=False, - base_units=( - self._reference_values if self._reference_values else None - ), + base_units=(self._reference_values if self._reference_values else None), ) self._snap_refs = self._reference_values.copy() return snap @@ -482,9 +460,7 @@ def _validate_forcefield(self, input_forcefield): else: # if there is only one ff for all molecule types ff_index = 0 - self._gmso_forcefields_dict[str(i)] = _force_field[ - ff_index - ].gmso_ff + self._gmso_forcefields_dict[str(i)] = _force_field[ff_index].gmso_ff def _assign_site_mol_type_idx(self): """Assign molecule type index to the gmso sites.""" @@ -563,24 +539,17 @@ def apply_forcefield( if not self._reference_values: epsilons = [ - s.atom_type.parameters["epsilon"] - for s in self.gmso_system.sites - ] - sigmas = [ - s.atom_type.parameters["sigma"] for s in self.gmso_system.sites + s.atom_type.parameters["epsilon"] for s in self.gmso_system.sites ] + sigmas = [s.atom_type.parameters["sigma"] for s in self.gmso_system.sites] masses = [s.mass for s in self.gmso_system.sites] energy_scale = np.max(epsilons) if self.auto_scale else 1.0 length_scale = np.max(sigmas) if self.auto_scale else 1.0 mass_scale = np.max(masses) if self.auto_scale else 1.0 - self._reference_values["energy"] = ( - energy_scale * epsilons[0].unit_array - ) - self._reference_values["length"] = ( - length_scale * sigmas[0].unit_array - ) + self._reference_values["energy"] = energy_scale * epsilons[0].unit_array + self._reference_values["length"] = length_scale * sigmas[0].unit_array self._reference_values["mass"] = mass_scale * masses[0].unit_array if remove_hydrogens: @@ -602,9 +571,7 @@ def visualize(self): if self.system: self.system.visualize().show() else: - raise ValueError( - "The initial configuraiton has not been created yet." - ) + raise ValueError("The initial configuraiton has not been created yet.") class Pack(System): @@ -660,8 +627,7 @@ def __init__( if not isinstance(density, u.array.unyt_quantity): self.density = density * u.Unit("g") / u.Unit("cm**3") warnings.warn( - "Units for density were not given, assuming " - "units of g/cm**3." + "Units for density were not given, assuming " "units of g/cm**3." ) else: self.density = density @@ -755,9 +721,7 @@ def _build_system(self): try: comp1 = self.all_molecules[next_idx] comp2 = self.all_molecules[next_idx + 1] - comp2.translate( - self.basis_vector * np.array([self.x, self.y, 0]) - ) + comp2.translate(self.basis_vector * np.array([self.x, self.y, 0])) # TODO: what if comp1 and comp2 have different names? unit_cell = mb.Compound( subcompounds=[comp1, comp2], name=comp1.name @@ -776,7 +740,5 @@ def _build_system(self): z_len = bounding_box.lengths[2] + 0.2 # Center the lattice in its box system.box = mb.box.Box(np.array([x_len, y_len, z_len])) - system.translate_to( - (system.box.Lx / 2, system.box.Ly / 2, system.box.Lz / 2) - ) + system.translate_to((system.box.Lx / 2, system.box.Ly / 2, system.box.Lz / 2)) return system diff --git a/flowermd/internal/ff_utils.py b/flowermd/internal/ff_utils.py index 932a0054..5d5ab8b0 100644 --- a/flowermd/internal/ff_utils.py +++ b/flowermd/internal/ff_utils.py @@ -20,9 +20,7 @@ def ff_xml_directory(): for dirpath, dirnames, filenames in os.walk(FF_DIR): for file in filenames: if file.endswith(".xml"): - ff_xml_directory[file.split(".xml")[0]] = os.path.join( - dirpath, file - ) + ff_xml_directory[file.split(".xml")[0]] = os.path.join(dirpath, file) return ff_xml_directory diff --git a/flowermd/internal/utils.py b/flowermd/internal/utils.py index a8f6bea1..3ebb44b7 100644 --- a/flowermd/internal/utils.py +++ b/flowermd/internal/utils.py @@ -77,9 +77,7 @@ def _is_float(num): raise ValueError("The reference value is not a number.") # if ref_value is an instance of unyt_quantity, check the dimension. - if isinstance(ref_value, u.unyt_quantity) and _is_valid_dimension( - ref_value.units - ): + if isinstance(ref_value, u.unyt_quantity) and _is_valid_dimension(ref_value.units): return ref_value # if ref_value is a string, check if it is a number and if it is, check if # the unit exists in unyt and has the correct dimension. diff --git a/flowermd/library/forcefields.py b/flowermd/library/forcefields.py index 3721ef8b..313218f7 100644 --- a/flowermd/library/forcefields.py +++ b/flowermd/library/forcefields.py @@ -61,9 +61,7 @@ class OPLS_AA_DIMETHYLETHER(BaseXMLForcefield): """OPLS All Atom for dimethyl ether molecule forcefield class.""" def __init__(self, forcefield_files=f"{FF_DIR}/dimethylether_opls.xml"): - super(OPLS_AA_DIMETHYLETHER, self).__init__( - forcefield_files=forcefield_files - ) + super(OPLS_AA_DIMETHYLETHER, self).__init__(forcefield_files=forcefield_files) self.description = ( "Based on flowermd.forcefields.OPLS_AA. " "Trimmed down to include only dimethyl ether parameters." @@ -408,9 +406,7 @@ def _load_file(file, **kwargs): pair_dict[pair_type]["U"] = table[:, 1] pair_dict[pair_type]["F"] = table[:, 2] if len(pair_r_min) != len(pair_r_max) != 1: - raise ValueError( - "All pair files must have the same r-range values" - ) + raise ValueError("All pair files must have the same r-range values") # Read bond files bond_dict = dict() if bonds: @@ -430,12 +426,9 @@ def _load_file(file, **kwargs): for angle_type in angles: table = _load_file(angles[angle_type], **kwargs) thetas = table[:, 0] - if thetas[0] != 0 or not np.allclose( - thetas[-1], np.pi, atol=1e-5 - ): + if thetas[0] != 0 or not np.allclose(thetas[-1], np.pi, atol=1e-5): raise ValueError( - "Angle values must be evenly spaced and " - "range from 0 to Pi." + "Angle values must be evenly spaced and " "range from 0 to Pi." ) angle_dict[angle_type] = dict() angle_dict[angle_type]["U"] = table[:, 1] @@ -446,9 +439,9 @@ def _load_file(file, **kwargs): for dih_type in dihedrals: table = _load_file(dihedrals[dih_type], **kwargs) thetas = table[:, 0] - if not np.allclose( - thetas[0], -np.pi, atol=1e-5 - ) or not np.allclose(thetas[-1], np.pi, atol=1e-5): + if not np.allclose(thetas[0], -np.pi, atol=1e-5) or not np.allclose( + thetas[-1], np.pi, atol=1e-5 + ): raise ValueError( "Dihedral angle values must be evenly spaced and " "range from -Pi to Pi." @@ -474,9 +467,7 @@ def _create_forcefield(self): nlist = hoomd.md.nlist.Cell( buffer=self.nlist_buffer, exclusions=self.exclusions ) - pair_table = hoomd.md.pair.Table( - nlist=nlist, default_r_cut=self.r_cut - ) + pair_table = hoomd.md.pair.Table(nlist=nlist, default_r_cut=self.r_cut) for pair_type in self.pairs: U = self.pairs[pair_type]["U"] F = self.pairs[pair_type]["F"] @@ -484,9 +475,7 @@ def _create_forcefield(self): raise ValueError( "The energy and force arrays are not the same size." ) - pair_table.params[tuple(pair_type)] = dict( - r_min=self.r_min, U=U, F=F - ) + pair_table.params[tuple(pair_type)] = dict(r_min=self.r_min, U=U, F=F) forces.append(pair_table) # Create bond forces if self.bonds: diff --git a/flowermd/library/polymers.py b/flowermd/library/polymers.py index 2ac9038f..ac99d311 100644 --- a/flowermd/library/polymers.py +++ b/flowermd/library/polymers.py @@ -258,9 +258,7 @@ def _build(self, length): bead_pair = "-".join([last_bead.name, next_bead.name]) bond_length = self.bond_lengths.get(bead_pair, None) if not bond_length: - bead_pair_rev = "-".join( - [next_bead.name, last_bead.name] - ) + bead_pair_rev = "-".join([next_bead.name, last_bead.name]) bond_length = self.bond_lengths.get(bead_pair_rev, None) if not bond_length: raise ValueError( @@ -319,12 +317,8 @@ def __init__(self, lengths, num_mols, lpar, bead_mass, bond_length): def _build(self, length): # Build up ellipsoid bead bead = mb.Compound(name="ellipsoid") - head = mb.Compound( - pos=(0, 0, self.lpar), name="A", mass=self.bead_mass / 4 - ) - tail = mb.Compound( - pos=(0, 0, -self.lpar), name="A", mass=self.bead_mass / 4 - ) + head = mb.Compound(pos=(0, 0, self.lpar), name="A", mass=self.bead_mass / 4) + tail = mb.Compound(pos=(0, 0, -self.lpar), name="A", mass=self.bead_mass / 4) head_mid = mb.Compound( pos=(0, 0, self.lpar / 2), name="B", mass=self.bead_mass / 4 ) diff --git a/flowermd/library/simulations/tensile.py b/flowermd/library/simulations/tensile.py index 7b2b37ea..c74d5b3c 100644 --- a/flowermd/library/simulations/tensile.py +++ b/flowermd/library/simulations/tensile.py @@ -66,16 +66,12 @@ def __init__( self.fix_right = hoomd.filter.Tags(right_tags.astype(np.uint32)) all_fixed = hoomd.filter.Union(self.fix_left, self.fix_right) # Set the group of particles to be integrated over - self.integrate_group = hoomd.filter.SetDifference( - hoomd.filter.All(), all_fixed - ) + self.integrate_group = hoomd.filter.SetDifference(hoomd.filter.All(), all_fixed) @property def strain(self): """The current strain of the simulation.""" - delta_L = ( - self.box_lengths_reduced[self._axis_index] - self.initial_length - ) + delta_L = self.box_lengths_reduced[self._axis_index] - self.initial_length return delta_L / self.initial_length def run_tensile(self, strain, n_steps, kT, tau_kt, period): diff --git a/flowermd/modules/surface_wetting/surface_wetting.py b/flowermd/modules/surface_wetting/surface_wetting.py index 3f35383e..e58f84fc 100644 --- a/flowermd/modules/surface_wetting/surface_wetting.py +++ b/flowermd/modules/surface_wetting/surface_wetting.py @@ -98,12 +98,11 @@ def run_droplet( """ # Shrink down to high density - if not isinstance( - shrink_density, u.array.unyt_quantity - ) and not isinstance(final_density, u.array.unyt_quantity): + if not isinstance(shrink_density, u.array.unyt_quantity) and not isinstance( + final_density, u.array.unyt_quantity + ): warnings.warn( - "Units for density were not given, assuming " - "units of g/cm**3." + "Units for density were not given, assuming " "units of g/cm**3." ) target_box_shrink = get_target_box_mass_density( density=shrink_density * (u.g / (u.cm**3)), @@ -123,10 +122,7 @@ def run_droplet( target_box_final = get_target_box_mass_density( density=final_density, mass=self.mass.to("g") ) - elif ( - shrink_density.units.dimensions - == number_density.units.dimensions - ): + elif shrink_density.units.dimensions == number_density.units.dimensions: raise ValueError( "For now, only mass density is supported " "in the surface wetting module." @@ -253,9 +249,7 @@ def _build_snapshot(self): ] self.drop_ptypes = self.drop_snapshot.particles.types - wetting_snapshot.particles.types = ( - self.surface_ptypes + self.drop_ptypes - ) + wetting_snapshot.particles.types = self.surface_ptypes + self.drop_ptypes wetting_snapshot.particles.typeid = np.concatenate( ( @@ -338,8 +332,7 @@ def _build_snapshot(self): self.surface_snapshot.dihedrals.N + self.drop_snapshot.dihedrals.N ) wetting_snapshot.dihedrals.types = ( - self.surface_snapshot.dihedrals.types - + self.drop_snapshot.dihedrals.types + self.surface_snapshot.dihedrals.types + self.drop_snapshot.dihedrals.types ) wetting_snapshot.dihedrals.typeid = np.concatenate( ( @@ -411,9 +404,7 @@ def _adjust_particle_positions(self, wetting_box): self.drop_snapshot.particles.position, axis=0 ) # shift drop particles z position to be at the top of surface - z_shift = ( - np.abs(min(drop_pos[:, 2]) - max(surface_pos[:, 2])) - self.gap - ) + z_shift = np.abs(min(drop_pos[:, 2]) - max(surface_pos[:, 2])) - self.gap drop_pos[:, 2] -= z_shift wetting_pos = np.concatenate((surface_pos, drop_pos), axis=0) return wetting_pos diff --git a/flowermd/modules/surface_wetting/utils.py b/flowermd/modules/surface_wetting/utils.py index c159530a..11a90dd9 100644 --- a/flowermd/modules/surface_wetting/utils.py +++ b/flowermd/modules/surface_wetting/utils.py @@ -102,9 +102,7 @@ def _combine_lj_forces( surface_ptype, drop_ptype, ) not in list(lj.params.keys()): - drop_epsilon = drop_lj.params[(drop_ptype, drop_ptype)][ - "epsilon" - ] + drop_epsilon = drop_lj.params[(drop_ptype, drop_ptype)]["epsilon"] surface_epsilon = surface_lj.params[ (surface_ptype[1:], surface_ptype[1:]) ]["epsilon"] diff --git a/flowermd/modules/welding/utils.py b/flowermd/modules/welding/utils.py index 029aa467..623f463e 100644 --- a/flowermd/modules/welding/utils.py +++ b/flowermd/modules/welding/utils.py @@ -50,13 +50,9 @@ def add_void_particles( ) # Set updated mass and charges init_mass = snapshot.particles.mass - snapshot.particles.mass = np.concatenate( - (init_mass, np.array([1])), axis=None - ) + snapshot.particles.mass = np.concatenate((init_mass, np.array([1])), axis=None) init_charges = snapshot.particles.charge - snapshot.particles.charge = np.concatenate( - (init_charges, np.array([0])), axis=None - ) + snapshot.particles.charge = np.concatenate((init_charges, np.array([0])), axis=None) # Updated LJ params lj = [i for i in forcefield if isinstance(i, hoomd.md.pair.LJ)][0] for ptype in snapshot.particles.types: diff --git a/flowermd/modules/welding/welding.py b/flowermd/modules/welding/welding.py index 3894b63d..f2518bc5 100644 --- a/flowermd/modules/welding/welding.py +++ b/flowermd/modules/welding/welding.py @@ -62,9 +62,7 @@ def _build(self): void_id = snap.particles.types.index("VOID") snap.particles.types.remove("VOID") keep_indices = np.where(snap.particles.typeid != void_id)[0] - snap.particles.position = snap.particles.position[ - keep_indices - ] + snap.particles.position = snap.particles.position[keep_indices] snap.particles.mass = snap.particles.mass[keep_indices] snap.particles.charge = snap.particles.charge[keep_indices] snap.particles.orientation = snap.particles.orientation[ @@ -89,18 +87,14 @@ def _build(self): # Set up snapshot.particles info: # Get set of new coordiantes, shifted along interface axis - shift = ( - snap_L.configuration.box[axis_index] + self.gap - self.wall_sigma - ) / 2 + shift = (snap_L.configuration.box[axis_index] + self.gap - self.wall_sigma) / 2 right_pos = np.copy(snap_R.particles.position) right_pos[:, axis_index] += shift left_pos = np.copy(snap_L.particles.position) left_pos[:, axis_index] -= shift pos = np.concatenate((left_pos, right_pos), axis=None) - mass = np.concatenate( - (snap_L.particles.mass, snap_R.particles.mass), axis=None - ) + mass = np.concatenate((snap_L.particles.mass, snap_R.particles.mass), axis=None) charges = np.concatenate( (snap_L.particles.charge, snap_R.particles.charge), axis=None ) @@ -116,9 +110,7 @@ def _build(self): # Set up bonds: bond_group_left = np.copy(snap_L.bonds.group) bond_group_right = np.copy(snap_R.bonds.group) + snap_R.particles.N - bond_group = np.concatenate( - (bond_group_left, bond_group_right), axis=None - ) + bond_group = np.concatenate((bond_group_left, bond_group_right), axis=None) bond_type_ids = np.concatenate( (snap_L.bonds.typeid, snap_R.bonds.typeid), axis=None ) @@ -129,9 +121,7 @@ def _build(self): # Set up angles: angle_group_left = np.copy(snap_L.angles.group) angle_group_right = np.copy(snap_R.angles.group) + snap_L.particles.N - angle_group = np.concatenate( - (angle_group_left, angle_group_right), axis=None - ) + angle_group = np.concatenate((angle_group_left, angle_group_right), axis=None) angle_type_ids = np.concatenate( (snap_L.angles.typeid, snap_R.angles.typeid), axis=None ) @@ -141,9 +131,7 @@ def _build(self): # Set up dihedrals: dihedral_group_left = np.copy(snap_L.dihedrals.group) - dihedral_group_right = ( - np.copy(snap_R.dihedrals.group) + snap_L.particles.N - ) + dihedral_group_right = np.copy(snap_R.dihedrals.group) + snap_L.particles.N dihedral_group = np.concatenate( (dihedral_group_left, dihedral_group_right), axis=None ) diff --git a/flowermd/tests/base/test_molecule.py b/flowermd/tests/base/test_molecule.py index 38c3b015..d0ef0d5a 100644 --- a/flowermd/tests/base/test_molecule.py +++ b/flowermd/tests/base/test_molecule.py @@ -43,9 +43,7 @@ def test_molecule_topology_benzene(self, benzene_mb): assert not any(molecule.topology_information["particle_charge"]) def test_validate_force_field_oplsaa(self, benzene_mb): - molecule = Molecule( - num_mols=2, force_field=OPLS_AA(), compound=benzene_mb - ) + molecule = Molecule(num_mols=2, force_field=OPLS_AA(), compound=benzene_mb) assert molecule.gmso_molecule.is_typed() assert set(molecule.topology_information["particle_types"]) == { "opls_145", @@ -66,23 +64,15 @@ def test_validate_force_field_xml_file_path(self, benzene_mb, benzene_xml): } assert any(molecule.topology_information["particle_charge"]) - def test_validate_force_field_hoomd_ff_aa( - self, benzene_mb, benzene_hoomd_ff - ): + def test_validate_force_field_hoomd_ff_aa(self, benzene_mb, benzene_hoomd_ff): hoomd_ff = benzene_hoomd_ff(include_hydrogen=True) - molecule = Molecule( - num_mols=2, force_field=hoomd_ff, compound=benzene_mb - ) + molecule = Molecule(num_mols=2, force_field=hoomd_ff, compound=benzene_mb) assert molecule.force_field == hoomd_ff assert not molecule.gmso_molecule.is_typed() - def test_validate_fore_field_hoomd_ff_ua( - self, benzene_mb, benzene_hoomd_ff - ): + def test_validate_fore_field_hoomd_ff_ua(self, benzene_mb, benzene_hoomd_ff): hoomd_ff = benzene_hoomd_ff(include_hydrogen=False) - molecule = Molecule( - num_mols=2, force_field=hoomd_ff, compound=benzene_mb - ) + molecule = Molecule(num_mols=2, force_field=hoomd_ff, compound=benzene_mb) assert molecule.force_field == hoomd_ff assert not molecule.gmso_molecule.is_typed() @@ -165,9 +155,7 @@ def test_coarse_grain_with_multiple_single_beads( assert molecule.molecules[0].n_particles == 4 assert molecule.topology_information["bond_types"] == {("A", "A")} - def test_coarse_grain_with_different_beads( - self, pps_smiles, benzene_smiles - ): + def test_coarse_grain_with_different_beads(self, pps_smiles, benzene_smiles): molecule = Molecule(num_mols=2, smiles=pps_smiles) molecule.coarse_grain(beads={"A": benzene_smiles, "B": "S"}) assert molecule.molecules[0].n_particles == 2 @@ -193,9 +181,7 @@ def test_polymer(self, dimethylether_smiles): assert polymer.n_particles == 23 assert polymer.n_bonds == 22 assert ("O", "C", "C") in polymer.topology_information["angle_types"] - assert ("O", "C", "C", "O") in polymer.topology_information[ - "dihedral_types" - ] + assert ("O", "C", "C", "O") in polymer.topology_information["dihedral_types"] def test_polymer_different_chain_lengths(self, dimethylether_smiles): polymer = Polymer( @@ -247,9 +233,7 @@ def test_copolymer_with_sequence(self, polyethylene, polyDME): ) assert copolymer.n_particles == 22 assert copolymer.random_sequence is False - assert ("C", "C", "C", "C") in copolymer.topology_information[ - "dihedral_types" - ] + assert ("C", "C", "C", "C") in copolymer.topology_information["dihedral_types"] def test_copolymer_with_sequence_different_chain_lengths( self, polyethylene, polyDME @@ -265,9 +249,7 @@ def test_copolymer_with_sequence_different_chain_lengths( assert copolymer.molecules[0].n_particles == 42 assert copolymer.molecules[1].n_particles == 62 - def test_copolymer_with_sequence_different_num_mol( - self, polyethylene, polyDME - ): + def test_copolymer_with_sequence_different_num_mol(self, polyethylene, polyDME): copolymer = CoPolymer( monomer_A=polyDME, monomer_B=polyethylene, @@ -290,9 +272,7 @@ def test_copolymer_random_sequence(self, polyethylene, polyDME): ) # sequence is BAA assert copolymer.n_particles == 22 - assert ("O", "C", "C", "O") in copolymer.topology_information[ - "dihedral_types" - ] + assert ("O", "C", "C", "O") in copolymer.topology_information["dihedral_types"] def test_align_z_axis(self, polyethylene): pe = polyethylene(num_mols=5, lengths=20) diff --git a/flowermd/tests/base/test_simulation.py b/flowermd/tests/base/test_simulation.py index a756f808..3c4b5271 100644 --- a/flowermd/tests/base/test_simulation.py +++ b/flowermd/tests/base/test_simulation.py @@ -27,9 +27,7 @@ def test_initialize_from_system(self, benzene_system): def test_initialize_from_system_separate_ff( self, benzene_cg_system, cg_single_bead_ff ): - sim = Simulation.from_system( - benzene_cg_system, forcefield=cg_single_bead_ff - ) + sim = Simulation.from_system(benzene_cg_system, forcefield=cg_single_bead_ff) sim.run_NVT(kT=0.1, tau_kt=10, n_steps=500) def test_initialize_from_system_missing_ff(self, benzene_cg_system): @@ -59,9 +57,7 @@ def test_reference_values(self, benzene_system): forcefield=benzene_system.hoomd_forcefield, reference_values=benzene_system.reference_values, ) - assert np.isclose( - float(sim.mass.value), benzene_system.mass.value, atol=1e-4 - ) + assert np.isclose(float(sim.mass.value), benzene_system.mass.value, atol=1e-4) assert np.allclose(benzene_system.box.lengths, sim.box_lengths.value) def test_set_ref_values(self, benzene_system): @@ -202,17 +198,13 @@ def test_update_volume_by_density_factor(self, benzene_system): period=1, final_box_lengths=target_box, ) - assert np.isclose( - sim.density.value, (init_density * 5).value, atol=1e-4 - ) + assert np.isclose(sim.density.value, (init_density * 5).value, atol=1e-4) def test_change_methods(self, benzene_system): sim = Simulation.from_system(benzene_system) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=0) assert isinstance(sim.method, hoomd.md.methods.ConstantVolume) - sim.run_NPT( - kT=1.0, tau_kt=0.01, tau_pressure=0.1, pressure=0.001, n_steps=0 - ) + sim.run_NPT(kT=1.0, tau_kt=0.01, tau_pressure=0.1, pressure=0.001, n_steps=0) assert isinstance(sim.method, hoomd.md.methods.ConstantPressure) def test_change_dt(self, benzene_system): @@ -404,9 +396,7 @@ def test_gsd_logger(self, benzene_system): def test_bad_ff(self, benzene_system): with pytest.raises(ValueError): - Simulation( - initial_state=benzene_system.hoomd_snapshot, forcefield="gaff" - ) + Simulation(initial_state=benzene_system.hoomd_snapshot, forcefield="gaff") with pytest.raises(ValueError): Simulation( initial_state=benzene_system.hoomd_snapshot, diff --git a/flowermd/tests/base/test_system.py b/flowermd/tests/base/test_system.py index cade5f96..5bb38de8 100644 --- a/flowermd/tests/base/test_system.py +++ b/flowermd/tests/base/test_system.py @@ -23,9 +23,7 @@ class TestSystem(BaseTest): def test_single_mol_type(self, benzene_molecule): benzene_mols = benzene_molecule(n_mols=3) system = Pack(molecules=[benzene_mols], density=0.8) - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(benzene_mols.molecules) assert system.gmso_system.is_typed() @@ -38,9 +36,7 @@ def test_multiple_mol_types(self, benzene_molecule, ethane_molecule): benzene_mol = benzene_molecule(n_mols=3) ethane_mol = ethane_molecule(n_mols=2) system = Pack(molecules=[benzene_mol, ethane_mol], density=0.8) - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) assert system.n_mol_types == 2 assert len(system.all_molecules) == len(benzene_mol.molecules) + len( ethane_mol.molecules @@ -69,9 +65,9 @@ def test_multiple_mol_types_different_ff( auto_scale=True, ) assert system.n_mol_types == 2 - assert len(system.all_molecules) == len( - dimethylether_mol.molecules - ) + len(pps_mol.molecules) + assert len(system.all_molecules) == len(dimethylether_mol.molecules) + len( + pps_mol.molecules + ) assert system.gmso_system.sites[0].group == "0" assert system.gmso_system.sites[-1].group == "1" assert system.gmso_system.is_typed() @@ -94,9 +90,7 @@ def test_multiple_mol_types_different_ff( def test_system_from_mol2_mol_parameterization(self, benzene_molecule_mol2): benzene_mol = benzene_molecule_mol2(n_mols=3) system = Pack(molecules=[benzene_mol], density=0.8) - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) assert system.gmso_system.is_typed() assert len(system.hoomd_forcefield) > 0 assert system.n_particles == system.hoomd_snapshot.particles.N @@ -110,9 +104,7 @@ def test_remove_hydrogen(self, benzene_molecule): remove_hydrogens=True, auto_scale=True, ) - assert not any( - [s.element.atomic_number == 1 for s in system.gmso_system.sites] - ) + assert not any([s.element.atomic_number == 1 for s in system.gmso_system.sites]) assert system.gmso_system.is_typed() assert len(system.hoomd_forcefield) > 0 assert list(system.hoomd_forcefield[0].params.keys()) == [ @@ -163,9 +155,7 @@ def test_remove_hydrogen_no_hydrogen(self, benzene_molecule): auto_scale=True, ) hydrogens = [ - site - for site in system.gmso_system.sites - if site.element.atomic_number == 1 + site for site in system.gmso_system.sites if site.element.atomic_number == 1 ] for h_site in hydrogens: system.gmso_system.remove_site(h_site) @@ -188,9 +178,7 @@ def test_add_mass_charges(self, benzene_molecule): assert site.charge == 0 snap = system.hoomd_snapshot - assert np.allclose( - sum(snap.particles.mass), 6 * (12.011 + 1.008), atol=1e-4 - ) + assert np.allclose(sum(snap.particles.mass), 6 * (12.011 + 1.008), atol=1e-4) assert sum(snap.particles.charge) == 0 def test_pack_box(self, benzene_molecule): @@ -213,13 +201,9 @@ def test_mass(self, pps_molecule): def test_ref_length(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) - assert np.allclose( - system.reference_length.to("angstrom").value, 3.5, atol=1e-3 - ) + assert np.allclose(system.reference_length.to("angstrom").value, 3.5, atol=1e-3) reduced_box = system.hoomd_snapshot.configuration.box[0:3] calc_box = reduced_box * system.reference_length.to("nm").value assert np.allclose(calc_box[0], system.box.Lx, atol=1e-2) @@ -229,9 +213,7 @@ def test_ref_length(self, polyethylene): def test_ref_mass(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) total_red_mass = sum(system.hoomd_snapshot.particles.mass) assert np.allclose( system.mass.value, @@ -242,9 +224,7 @@ def test_ref_mass(self, polyethylene): def test_ref_energy(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) assert np.allclose( system.reference_energy.to("kcal/mol").value, 0.066, atol=1e-3 ) @@ -255,9 +235,7 @@ def test_rebuild_snapshot(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) assert system._snap_refs == system.reference_values assert system._ff_refs == system.reference_values init_snap = system.hoomd_snapshot @@ -277,9 +255,7 @@ def test_ref_values_no_autoscale(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_length = 1 * u.angstrom system.reference_energy = 1 * u.kcal / u.mol system.reference_mass = 1 * u.amu @@ -289,16 +265,12 @@ def test_ref_values_no_autoscale(self, polyethylene): ) assert dict(system.hoomd_forcefield[3].params)["opls_135", "opls_135"][ "epsilon" - ] == system.gmso_system.sites[0].atom_type.parameters["epsilon"].to( - "kcal/mol" - ) + ] == system.gmso_system.sites[0].atom_type.parameters["epsilon"].to("kcal/mol") def test_set_ref_values(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=1) system = Pack(molecules=[polyethylene], density=1.0) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -315,9 +287,7 @@ def test_set_ref_values_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) ref_value_dict = { "length": "1 angstrom", "energy": "3 kcal/mol", @@ -334,9 +304,7 @@ def test_set_ref_values_missing_key(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -350,9 +318,7 @@ def test_set_ref_values_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -367,9 +333,7 @@ def test_set_ref_values_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -384,9 +348,7 @@ def test_set_ref_length(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_length = 1 * u.angstrom assert system.reference_length == 1 * u.angstrom @@ -396,9 +358,7 @@ def test_set_ref_length_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_length = 1.0 @@ -408,9 +368,7 @@ def test_ref_length_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_length = "1 angstrom" assert system.reference_length == 1 * u.angstrom @@ -420,9 +378,7 @@ def test_ref_length_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0" @@ -432,9 +388,7 @@ def test_ref_length_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0 invalid_unit" @@ -444,9 +398,7 @@ def test_ref_length_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_length = 1.0 * u.g @@ -456,9 +408,7 @@ def test_ref_length_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0 g" @@ -468,9 +418,7 @@ def test_ref_length_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) system.reference_length = 1 * u.angstrom assert system.reference_length == 1 * u.angstrom @@ -480,9 +428,7 @@ def test_set_ref_energy(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_energy = 1 * u.kcal / u.mol assert system.reference_energy == 1 * u.kcal / u.mol @@ -492,9 +438,7 @@ def test_set_ref_energy_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 @@ -504,9 +448,7 @@ def test_ref_energy_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_energy = "1 kJ" assert system.reference_energy == 1 * u.kJ @@ -516,9 +458,7 @@ def test_ref_energy_string_comb_units(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_energy = "1 kcal/mol" assert system.reference_energy == 1 * u.kcal / u.mol @@ -528,9 +468,7 @@ def test_ref_energy_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0" @@ -540,9 +478,7 @@ def test_ref_energy_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0 invalid_unit" @@ -552,9 +488,7 @@ def test_ref_energy_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 * u.g @@ -564,9 +498,7 @@ def test_ref_energy_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0 m" @@ -576,9 +508,7 @@ def test_set_ref_energy_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) system.reference_energy = 1 * u.kcal / u.mol assert system.reference_energy == 1 * u.kcal / u.mol @@ -588,9 +518,7 @@ def test_set_ref_mass(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_mass = 1.0 * u.amu assert system.reference_mass == 1.0 * u.amu @@ -601,9 +529,7 @@ def test_set_ref_mass_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_mass = 1.0 @@ -613,9 +539,7 @@ def test_ref_mass_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.reference_mass = "1 g" assert system.reference_mass == 1.0 * u.g @@ -625,9 +549,7 @@ def test_ref_mass_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0" @@ -637,9 +559,7 @@ def test_ref_mass_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0 invalid_unit" @@ -649,9 +569,7 @@ def test_ref_mass_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 * u.m @@ -661,9 +579,7 @@ def test_ref_mass_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0 m" @@ -673,9 +589,7 @@ def test_set_ref_mass_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) system.reference_mass = 1.0 * u.amu assert system.reference_mass == 1.0 * u.amu @@ -687,9 +601,7 @@ def test_apply_forcefield_no_forcefield(self, polyethylene): density=1.0, ) with pytest.raises(ForceFieldError): - system.apply_forcefield( - r_cut=2.5, force_field=None, auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=None, auto_scale=False) def test_apply_forcefield_no_forcefield_w_mol_ff(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=3, force_field=OPLS_AA()) @@ -709,9 +621,7 @@ def test_apply_forcefield_w_mol_ff(self, benzene_molecule): density=1.0, ) with pytest.raises(ForceFieldError): - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) def test_validate_forcefield_invalid_ff_type(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=1) @@ -762,9 +672,7 @@ def test_forcefield_list_hoomd_ff(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) hoomd_ff = system.hoomd_forcefield @@ -783,9 +691,7 @@ def test_lattice_polymer(self, polyethylene): y=1, n=4, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(polyethylene.molecules) @@ -795,9 +701,7 @@ def test_lattice_polymer(self, polyethylene): for mol_class in system._molecules: for mol in mol_class.molecules: backbone = get_backbone_vector(mol.xyz) - assert np.allclose( - np.abs(backbone), np.array([0, 0, 1]), atol=1e-1 - ) + assert np.allclose(np.abs(backbone), np.array([0, 0, 1]), atol=1e-1) def test_lattice_molecule(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=32) @@ -807,9 +711,7 @@ def test_lattice_molecule(self, benzene_molecule): y=1, n=4, ) - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(benzene_mol.molecules) assert len(system.hoomd_forcefield) > 0 @@ -852,9 +754,7 @@ def test_to_gsd(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.to_gsd("test.gsd") assert os.path.isfile(os.path.join(os.getcwd(), "test.gsd")) os.remove(os.path.join(os.getcwd(), "test.gsd")) @@ -875,9 +775,7 @@ def test_pickle_ff(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) system.pickle_forcefield("forcefield.pickle") assert os.path.isfile(os.path.join(os.getcwd(), "forcefield.pickle")) os.remove(os.path.join(os.getcwd(), "forcefield.pickle")) diff --git a/flowermd/tests/base_test.py b/flowermd/tests/base_test.py index b9637c3c..3b8b7867 100644 --- a/flowermd/tests/base_test.py +++ b/flowermd/tests/base_test.py @@ -103,9 +103,7 @@ def _hoomd_ff(include_hydrogen, invalid_pair=False): ) bonds = self.benzene_hoomd_bond(include_hydrogen=include_hydrogen) angles = self.benzene_hoomd_angle(include_hydrogen=include_hydrogen) - dihedrals = self.benzene_hoomd_dihedral( - include_hydrogen=include_hydrogen - ) + dihedrals = self.benzene_hoomd_dihedral(include_hydrogen=include_hydrogen) return [pairs, bonds, angles, dihedrals] return _hoomd_ff @@ -147,9 +145,7 @@ def _pps_molecule(n_mols): @pytest.fixture() def dimethylether_molecule(self, dimethylether_smiles): def _dimethylether_molecule(n_mols): - dimethylether = Molecule( - num_mols=n_mols, smiles=dimethylether_smiles - ) + dimethylether = Molecule(num_mols=n_mols, smiles=dimethylether_smiles) return dimethylether return _dimethylether_molecule @@ -169,7 +165,7 @@ def __init__(self, lengths, num_mols, **kwargs): bond_indices=bond_indices, bond_length=bond_length, bond_orientation=bond_orientation, - **kwargs + **kwargs, ) return _PolyEthylene @@ -189,7 +185,7 @@ def __init__(self, lengths, num_mols, **kwargs): bond_indices=bond_indices, bond_length=bond_length, bond_orientation=bond_orientation, - **kwargs + **kwargs, ) return _PPS @@ -209,7 +205,7 @@ def __init__(self, lengths, num_mols, **kwargs): bond_indices=bond_indices, bond_length=bond_length, bond_orientation=bond_orientation, - **kwargs + **kwargs, ) return _PolyDME @@ -221,9 +217,7 @@ def benzene_system(self, benzene_mb): molecules=[benzene], density=0.2 * u.g / u.cm**3, ) - system.apply_forcefield( - r_cut=2.5, force_field=OPLS_AA(), auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) return system @pytest.fixture() diff --git a/flowermd/tests/library/test_forcefields.py b/flowermd/tests/library/test_forcefields.py index 30ae0920..fb927109 100644 --- a/flowermd/tests/library/test_forcefields.py +++ b/flowermd/tests/library/test_forcefields.py @@ -117,15 +117,9 @@ def test_ellipsoid_ff(self): angle_theta0=1.57, ) assert len(ellipsoid_ff.hoomd_forces) == 3 - assert isinstance( - ellipsoid_ff.hoomd_forces[-1], hoomd.md.pair.aniso.GayBerne - ) - assert ("R", "R") in list( - dict(ellipsoid_ff.hoomd_forces[-1].params).keys() - ) - assert ("B", "R") in list( - dict(ellipsoid_ff.hoomd_forces[-1].params).keys() - ) + assert isinstance(ellipsoid_ff.hoomd_forces[-1], hoomd.md.pair.aniso.GayBerne) + assert ("R", "R") in list(dict(ellipsoid_ff.hoomd_forces[-1].params).keys()) + assert ("B", "R") in list(dict(ellipsoid_ff.hoomd_forces[-1].params).keys()) assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["epsilon"] == 0.0 assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["lperp"] == 0.0 assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["lpar"] == 0.0 @@ -155,9 +149,7 @@ def test_from_txt_file(self): def test_from_csv_file(self): pair_file = os.path.join(ASSETS_DIR, "lj_pair_table.csv") - ff = TableForcefield.from_files( - pairs={("A", "A"): pair_file}, delimiter="," - ) + ff = TableForcefield.from_files(pairs={("A", "A"): pair_file}, delimiter=",") pair_data = np.genfromtxt(pair_file, delimiter=",") assert ff.r_min == pair_data[:, 0][0] assert ff.r_cut == pair_data[:, 0][-1] diff --git a/flowermd/tests/library/test_tensile.py b/flowermd/tests/library/test_tensile.py index 161dfc38..7cedc9ec 100644 --- a/flowermd/tests/library/test_tensile.py +++ b/flowermd/tests/library/test_tensile.py @@ -14,9 +14,7 @@ def test_tensile(self): y=1.2, n=4, ) - system.apply_forcefield( - r_cut=2.5, force_field=[OPLS_AA_PPS()], auto_scale=True - ) + system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA_PPS()], auto_scale=True) tensile_sim = Tensile( initial_state=system.hoomd_snapshot, diff --git a/flowermd/tests/modules/test_surface_wetting.py b/flowermd/tests/modules/test_surface_wetting.py index b2dc908b..345beda9 100644 --- a/flowermd/tests/modules/test_surface_wetting.py +++ b/flowermd/tests/modules/test_surface_wetting.py @@ -43,9 +43,7 @@ def test_droplet_sim(self, polyethylene): final_density=0.05 * u.g / u.cm**3, tau_kt=drop_sim.dt * 100, ) - assert np.isclose( - drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 - ) + assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) def test_droplet_sim_no_units(self, polyethylene): drop_mol = polyethylene(num_mols=10, lengths=5) @@ -76,9 +74,7 @@ def test_droplet_sim_no_units(self, polyethylene): final_density=0.05, tau_kt=drop_sim.dt * 100, ) - assert np.isclose( - drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 - ) + assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) def test_droplet_sim_bad_units(self, polyethylene): drop_mol = polyethylene(num_mols=10, lengths=5) @@ -109,9 +105,7 @@ def test_droplet_sim_bad_units(self, polyethylene): final_density=0.05 * (u.cm**-3), tau_kt=drop_sim.dt * 100, ) - assert np.isclose( - drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 - ) + assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) class TestInterfaceBuilder(BaseTest): @@ -160,8 +154,7 @@ def test_interface_builder( + drop_snapshot.particles.types ) assert np.isclose( - interface.hoomd_snapshot.configuration.box[2] - * drop_refs["length"].value, + interface.hoomd_snapshot.configuration.box[2] * drop_refs["length"].value, 15, atol=1e-2, ) @@ -189,9 +182,7 @@ def test_interface_builder( class TestWettingSimulation(BaseTest): - def test_wetting_sim( - self, surface_wetting_init_snapshot, surface_wetting_init_ff - ): + def test_wetting_sim(self, surface_wetting_init_snapshot, surface_wetting_init_ff): # load surface wetting init snapshot snapshot = gsd.hoomd.open(surface_wetting_init_snapshot)[0] # load ff from pickle diff --git a/flowermd/tests/modules/test_weld.py b/flowermd/tests/modules/test_weld.py index 09ed0e42..06f87ec7 100644 --- a/flowermd/tests/modules/test_weld.py +++ b/flowermd/tests/modules/test_weld.py @@ -90,9 +90,7 @@ def test_slab_sim_xaxis(self, polyethylene_system): forcefield=polyethylene_system.hoomd_forcefield, log_write_freq=2000, ) - assert any( - [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] - ) + assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_slab_sim_yaxis(self, polyethylene_system): @@ -102,9 +100,7 @@ def test_slab_sim_yaxis(self, polyethylene_system): interface_axis=(0, 1, 0), log_write_freq=2000, ) - assert any( - [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] - ) + assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_slab_sim_zaxis(self, polyethylene_system): @@ -114,9 +110,7 @@ def test_slab_sim_zaxis(self, polyethylene_system): interface_axis=(0, 0, 1), log_write_freq=2000, ) - assert any( - [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] - ) + assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_weld_sim(self, polyethylene_system): @@ -128,9 +122,7 @@ def test_weld_sim(self, polyethylene_system): sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) sim.save_restart_gsd() # Create interface system from slab restart.gsd file - interface = Interface( - gsd_files="restart.gsd", interface_axis="x", gap=0.1 - ) + interface = Interface(gsd_files="restart.gsd", interface_axis="x", gap=0.1) sim = WeldSimulation( initial_state=interface.hoomd_snapshot, forcefield=polyethylene_system.hoomd_forcefield, diff --git a/flowermd/tests/utils/test_actions.py b/flowermd/tests/utils/test_actions.py index 3955bfc3..7cb2f90d 100644 --- a/flowermd/tests/utils/test_actions.py +++ b/flowermd/tests/utils/test_actions.py @@ -11,24 +11,18 @@ class TestActions(BaseTest): def test_scale_epsilon(self, benzene_simulation): sim = benzene_simulation epsilon_scale = ScaleEpsilon(sim=sim, scale_factor=0.5) - energy_operation = hoomd.update.CustomUpdater( - action=epsilon_scale, trigger=10 - ) + energy_operation = hoomd.update.CustomUpdater(action=epsilon_scale, trigger=10) sim.operations.updaters.append(energy_operation) old_lj_force = copy.deepcopy(sim._lj_force().params) sim.run_NVT(n_steps=10, kT=1.0, tau_kt=1.0) new_lj_force = sim._lj_force().params for k in old_lj_force.keys(): - assert ( - new_lj_force[k]["epsilon"] == old_lj_force[k]["epsilon"] + 0.5 - ) + assert new_lj_force[k]["epsilon"] == old_lj_force[k]["epsilon"] + 0.5 def test_scale_sigma(self, benzene_simulation): sim = benzene_simulation sigma_scale = ScaleSigma(sim=sim, scale_factor=0.5) - energy_operation = hoomd.update.CustomUpdater( - action=sigma_scale, trigger=10 - ) + energy_operation = hoomd.update.CustomUpdater(action=sigma_scale, trigger=10) sim.operations.updaters.append(energy_operation) old_lj_force = copy.deepcopy(sim._lj_force().params) sim.run_NVT(n_steps=10, kT=1.0, tau_kt=1.0) diff --git a/flowermd/tests/utils/test_rigid_body.py b/flowermd/tests/utils/test_rigid_body.py index 08ed20d5..9c5f5c52 100644 --- a/flowermd/tests/utils/test_rigid_body.py +++ b/flowermd/tests/utils/test_rigid_body.py @@ -29,8 +29,7 @@ def test_ellipsoid_create_rigid_body(self): ) assert rigid_frame.particles.N == 8 + ellipsoid_chain.n_particles assert ( - rigid_frame.particles.types - == ["R"] + system.hoomd_snapshot.particles.types + rigid_frame.particles.types == ["R"] + system.hoomd_snapshot.particles.types ) assert rigid_frame.particles.mass[0] == 100 assert np.all( diff --git a/flowermd/tests/utils/test_utils.py b/flowermd/tests/utils/test_utils.py index 505eb0d3..06283edd 100644 --- a/flowermd/tests/utils/test_utils.py +++ b/flowermd/tests/utils/test_utils.py @@ -30,9 +30,7 @@ def test_validate_ref_value(self): assert validate_ref_value("1.0 kcal/mol", u.dimensions.energy) == ( 1.0 * u.kcal / u.mol ) - assert validate_ref_value("1.0 amu", u.dimensions.mass) == 1.0 * u.Unit( - "amu" - ) + assert validate_ref_value("1.0 amu", u.dimensions.mass) == 1.0 * u.Unit("amu") with pytest.raises(ReferenceUnitError): validate_ref_value("1.0 g", u.dimensions.energy) @@ -81,9 +79,7 @@ def test_target_box_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - target_box = get_target_box_number_density( - density=density, n_beads=n_beads - ) + target_box = get_target_box_number_density(density=density, n_beads=n_beads) L = target_box[0].value assert np.allclose(L**3, 100, atol=1e-8) @@ -91,9 +87,7 @@ def test_target_box_one_constraint_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - cubic_box = get_target_box_number_density( - density=density, n_beads=n_beads - ) + cubic_box = get_target_box_number_density(density=density, n_beads=n_beads) tetragonal_box = get_target_box_number_density( density=density, n_beads=n_beads, @@ -106,9 +100,7 @@ def test_target_box_two_constraint_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - cubic_box = get_target_box_number_density( - density=density, n_beads=n_beads - ) + cubic_box = get_target_box_number_density(density=density, n_beads=n_beads) ortho_box = get_target_box_number_density( density=density, n_beads=n_beads, @@ -116,9 +108,7 @@ def test_target_box_two_constraint_number_density(self): y_constraint=cubic_box[0] / 2, ) assert cubic_box[0] / 2 == ortho_box[0] == ortho_box[1] != ortho_box[2] - assert np.allclose( - ortho_box[2].value, cubic_box[0].value * 4, atol=1e-5 - ) + assert np.allclose(ortho_box[2].value, cubic_box[0].value * 4, atol=1e-5) def test_calculate_box_length_bad_args(self): mass_density = 1 * u.g / (u.cm**3) @@ -132,16 +122,12 @@ def test_calculate_box_length_fixed_l_1d(self): mass = u.unyt_quantity(6.0, u.g) density = u.unyt_quantity(0.5, u.g / u.cm**3) fixed_L = u.unyt_quantity(3.0, u.cm) - box_length = _calculate_box_length( - mass=mass, density=density, fixed_L=fixed_L - ) + box_length = _calculate_box_length(mass=mass, density=density, fixed_L=fixed_L) assert box_length == 2.0 * u.cm def test_calculate_box_length_fixed_l_2d(self): mass = u.unyt_quantity(12.0, u.g) density = u.unyt_quantity(0.5, u.g / u.cm**3) fixed_L = u.unyt_array([3.0, 2.0], u.cm) - box_length = _calculate_box_length( - mass=mass, density=density, fixed_L=fixed_L - ) + box_length = _calculate_box_length(mass=mass, density=density, fixed_L=fixed_L) assert box_length == 4.0 * u.cm diff --git a/flowermd/utils/rigid_body.py b/flowermd/utils/rigid_body.py index 4484b1e6..fba8bb36 100644 --- a/flowermd/utils/rigid_body.py +++ b/flowermd/utils/rigid_body.py @@ -70,9 +70,7 @@ def create_rigid_body( """ # find typeid sequence of the constituent particles types in a rigid bead p_types = np.array(snapshot.particles.types) - constituent_type_ids = np.where( - p_types[:, None] == bead_constituents_types - )[0] + constituent_type_ids = np.where(p_types[:, None] == bead_constituents_types)[0] # find indices that matches the constituent particle types typeids = snapshot.particles.typeid @@ -88,9 +86,7 @@ def create_rigid_body( n_rigid = rigid_const_idx.shape[0] # number of rigid bodies # calculate center of mass and its position for each rigid body - com_mass, com_position, com_moi = _get_com_mass_pos_moi( - snapshot, rigid_const_idx - ) + com_mass, com_position, com_moi = _get_com_mass_pos_moi(snapshot, rigid_const_idx) rigid_frame = gsd.hoomd.Frame() rigid_frame.particles.types = [bead_name] + snapshot.particles.types @@ -98,9 +94,7 @@ def create_rigid_body( rigid_frame.particles.typeid = np.concatenate( (([0] * n_rigid), snapshot.particles.typeid + 1) ) - rigid_frame.particles.mass = np.concatenate( - (com_mass, snapshot.particles.mass) - ) + rigid_frame.particles.mass = np.concatenate((com_mass, snapshot.particles.mass)) rigid_frame.particles.position = np.concatenate( (com_position, snapshot.particles.position) ) @@ -146,9 +140,7 @@ def create_rigid_body( # find local coordinates of the particles in the first rigid body # only need to find the local coordinates for the first rigid body - local_coords = ( - snapshot.particles.position[rigid_const_idx[0]] - com_position[0] - ) + local_coords = snapshot.particles.position[rigid_const_idx[0]] - com_position[0] rigid_constrain = hoomd.md.constrain.Rigid() rigid_constrain.body["R"] = { diff --git a/flowermd/utils/utils.py b/flowermd/utils/utils.py index afb08949..e2c973ab 100644 --- a/flowermd/utils/utils.py +++ b/flowermd/utils/utils.py @@ -88,9 +88,7 @@ def get_target_box_number_density( else: constraints = np.array([x_constraint, y_constraint, z_constraint]) fixed_L = constraints[np.not_equal(constraints, None).nonzero()] - L = _calculate_box_length( - density=density, n_beads=n_beads, fixed_L=fixed_L - ) + L = _calculate_box_length(density=density, n_beads=n_beads, fixed_L=fixed_L) constraints[np.equal(constraints, None).nonzero()] = L Lx, Ly, Lz = constraints return np.array([Lx, Ly, Lz]) * u.Unit(Lx.units) diff --git a/setup.py b/setup.py index cd3f223a..1b717a86 100644 --- a/setup.py +++ b/setup.py @@ -6,9 +6,7 @@ # Package meta-data. NAME = "flowermd" -DESCRIPTION = ( - "Package making it easier to build and simulate polymers in Hoomd-Blue." -) +DESCRIPTION = "Package making it easier to build and simulate polymers in Hoomd-Blue." URL = "https://github.com/cmelab/flowerMD" EMAIL = "chrisjones4@u.boisestate.edu" AUTHOR = "CME Lab" @@ -58,9 +56,7 @@ def run(self): pass self.status("Building Source and Wheel (universal) distribution…") - os.system( - "{0} setup.py sdist bdist_wheel --universal".format(sys.executable) - ) + os.system("{0} setup.py sdist bdist_wheel --universal".format(sys.executable)) self.status("Uploading the package to PyPi via Twine…") os.system("twine upload dist/*") From a2b95b3b5cd39d6b3c0c5500c35a82ceea12a32e Mon Sep 17 00:00:00 2001 From: chrisjonesbsu Date: Wed, 8 May 2024 14:53:21 -0600 Subject: [PATCH 2/3] remove import * --- flowermd/__init__.py | 1 + flowermd/assets/__init__.py | 1 + flowermd/base/__init__.py | 1 + flowermd/internal/__init__.py | 1 + flowermd/library/__init__.py | 1 + flowermd/modules/surface_wetting/__init__.py | 1 + flowermd/modules/welding/__init__.py | 1 + flowermd/tests/__init__.py | 1 + flowermd/utils/__init__.py | 11 ++++++++++- 9 files changed, 18 insertions(+), 1 deletion(-) diff --git a/flowermd/__init__.py b/flowermd/__init__.py index e0995a2c..f8c66533 100644 --- a/flowermd/__init__.py +++ b/flowermd/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """flowerMD package.""" from .base import ( diff --git a/flowermd/assets/__init__.py b/flowermd/assets/__init__.py index 7a63aa3f..dd30075b 100644 --- a/flowermd/assets/__init__.py +++ b/flowermd/assets/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """Paths to the assets used by flowerMD.""" from .forcefields import FF_DIR diff --git a/flowermd/base/__init__.py b/flowermd/base/__init__.py index 4376da88..14c369e0 100644 --- a/flowermd/base/__init__.py +++ b/flowermd/base/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """Base classes for flowerMD.""" from .forcefield import BaseHOOMDForcefield, BaseXMLForcefield diff --git a/flowermd/internal/__init__.py b/flowermd/internal/__init__.py index d449c3b3..a06e7713 100644 --- a/flowermd/internal/__init__.py +++ b/flowermd/internal/__init__.py @@ -1,2 +1,3 @@ +# ruff: noqa: F401 from .ff_utils import xml_to_gmso_ff from .utils import check_return_iterable, validate_ref_value diff --git a/flowermd/library/__init__.py b/flowermd/library/__init__.py index 64ea7c0d..d57a1bc2 100644 --- a/flowermd/library/__init__.py +++ b/flowermd/library/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """Library of predefined molecules, recipes and forcefields.""" from .forcefields import ( diff --git a/flowermd/modules/surface_wetting/__init__.py b/flowermd/modules/surface_wetting/__init__.py index b5a838ef..2c194ea0 100644 --- a/flowermd/modules/surface_wetting/__init__.py +++ b/flowermd/modules/surface_wetting/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """Surface wetting module for FlowerMD.""" from .surface_wetting import ( diff --git a/flowermd/modules/welding/__init__.py b/flowermd/modules/welding/__init__.py index d927c0a6..4768054a 100644 --- a/flowermd/modules/welding/__init__.py +++ b/flowermd/modules/welding/__init__.py @@ -1,3 +1,4 @@ +# ruff: noqa: F401 """Welding module for FlowerMD.""" from .utils import add_void_particles diff --git a/flowermd/tests/__init__.py b/flowermd/tests/__init__.py index 9eb4c5fb..2bb6e9cb 100644 --- a/flowermd/tests/__init__.py +++ b/flowermd/tests/__init__.py @@ -1 +1,2 @@ +# ruff: noqa: F401 from .base_test import BaseTest diff --git a/flowermd/utils/__init__.py b/flowermd/utils/__init__.py index f3e56cb6..5915613c 100644 --- a/flowermd/utils/__init__.py +++ b/flowermd/utils/__init__.py @@ -1,4 +1,13 @@ -from .actions import * +# ruff: noqa: F401 +"""Helpful utility functions for use with flowerMD.""" + +from .actions import ( + PullParticles, + ScaleEpsilon, + ScaleSigma, + StdOutLogger, + UpdateWalls, +) from .base_types import HOOMDThermostats from .rigid_body import create_rigid_body from .utils import ( From 3bb186f6ee37eb3f072d442c4f9d6f2f838a0efa Mon Sep 17 00:00:00 2001 From: chrisjonesBSU Date: Wed, 15 May 2024 20:32:45 -0600 Subject: [PATCH 3/3] change max line length back to 80 --- .pre-commit-config.yaml | 3 +- docs/source/conf.py | 4 +- flowermd/base/forcefield.py | 8 +- flowermd/base/molecule.py | 21 +- flowermd/base/simulation.py | 47 +++- flowermd/base/system.py | 82 +++++-- flowermd/internal/ff_utils.py | 4 +- flowermd/internal/utils.py | 4 +- flowermd/library/forcefields.py | 29 ++- flowermd/library/polymers.py | 12 +- flowermd/library/simulations/tensile.py | 8 +- .../surface_wetting/surface_wetting.py | 25 ++- flowermd/modules/surface_wetting/utils.py | 4 +- flowermd/modules/welding/utils.py | 8 +- flowermd/modules/welding/welding.py | 24 +- flowermd/tests/base/test_molecule.py | 40 +++- flowermd/tests/base/test_simulation.py | 20 +- flowermd/tests/base/test_system.py | 210 +++++++++++++----- flowermd/tests/base_test.py | 12 +- flowermd/tests/library/test_forcefields.py | 16 +- flowermd/tests/library/test_tensile.py | 4 +- .../tests/modules/test_surface_wetting.py | 19 +- flowermd/tests/modules/test_weld.py | 16 +- flowermd/tests/utils/test_actions.py | 12 +- flowermd/tests/utils/test_rigid_body.py | 3 +- flowermd/tests/utils/test_utils.py | 28 ++- flowermd/utils/rigid_body.py | 16 +- flowermd/utils/utils.py | 4 +- setup.py | 8 +- 29 files changed, 511 insertions(+), 180 deletions(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index c6c5daf3..6461d9c5 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -12,8 +12,9 @@ repos: rev: v0.4.2 # Ruff version hooks: - id: ruff - args: [--line-length=80, --fix, --extend-ignore=E203] + args: [--fix, --extend-ignore=E203] - id: ruff-format + args: [ --line-length=80 ] - repo: https://github.com/pre-commit/pre-commit-hooks rev: v4.6.0 hooks: diff --git a/docs/source/conf.py b/docs/source/conf.py index a6ccd74a..18e21b53 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -10,7 +10,9 @@ import sys project = "flowerMD" -copyright = "2023, Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" +copyright = ( + "2023, Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" +) author = "Chris Jones, Marjan Albooyeh, Rainier Barrett, Eric Jankowski" sys.path.insert(0, os.path.abspath("../..")) diff --git a/flowermd/base/forcefield.py b/flowermd/base/forcefield.py index 62dff440..532669d4 100644 --- a/flowermd/base/forcefield.py +++ b/flowermd/base/forcefield.py @@ -11,7 +11,9 @@ def __init__(self, forcefield_files=None, name=None): super(BaseXMLForcefield, self).__init__( forcefield_files=forcefield_files, name=name ) - self.gmso_ff = ffutils.FoyerFFs().load(forcefield_files or name).to_gmso_ff() + self.gmso_ff = ( + ffutils.FoyerFFs().load(forcefield_files or name).to_gmso_ff() + ) class BaseHOOMDForcefield: @@ -20,6 +22,8 @@ class BaseHOOMDForcefield: def __init__(self, hoomd_forces): self.hoomd_forces = hoomd_forces if hoomd_forces is None: - raise NotImplementedError("`hoomd_forces` must be defined in the subclass.") + raise NotImplementedError( + "`hoomd_forces` must be defined in the subclass." + ) if not isinstance(hoomd_forces, list): raise TypeError("`hoomd_forces` must be a list.") diff --git a/flowermd/base/molecule.py b/flowermd/base/molecule.py index a37d2a3c..298f5f1d 100644 --- a/flowermd/base/molecule.py +++ b/flowermd/base/molecule.py @@ -204,7 +204,8 @@ def _load(self): return mb.load(self.smiles, smiles=True) else: raise MoleculeLoadError( - msg=f"Unable to load the molecule from smiles " f"{self.smiles}." + msg=f"Unable to load the molecule from smiles " + f"{self.smiles}." ) def _align_backbones_z_axis(self, heavy_atoms_only=False): @@ -213,7 +214,11 @@ def _align_backbones_z_axis(self, heavy_atoms_only=False): if heavy_atoms_only: try: positions = np.array( - [p.xyz[0] for p in mol.particles() if p.element.symbol != "H"] + [ + p.xyz[0] + for p in mol.particles() + if p.element.symbol != "H" + ] ) except AttributeError: positions = mol.xyz @@ -272,7 +277,9 @@ def _identify_particle_information(self, gmso_molecule): ): self.hydrogen_types.append(p_name) self.particle_typeid.append(self.particle_types.index(p_name)) - self.particle_charge.append(site.charge.to_value() if site.charge else 0) + self.particle_charge.append( + site.charge.to_value() if site.charge else 0 + ) def _identify_pairs(self, particle_types): """Identify all unique particle pairs from the particle types. @@ -283,7 +290,9 @@ def _identify_pairs(self, particle_types): List of all particle types. """ - self.pairs = set(itertools.combinations_with_replacement(particle_types, 2)) + self.pairs = set( + itertools.combinations_with_replacement(particle_types, 2) + ) def _identify_bond_types(self, gmso_molecule): """Identify all unique bond types from the GMSO topology. @@ -426,7 +435,9 @@ def _validate_force_field(self): # Update topology information from typed gmso after applying ff. self._identify_topology_information(self.gmso_molecule) elif isinstance(self.force_field, BaseHOOMDForcefield): - _validate_hoomd_ff(self.force_field.hoomd_forces, self.topology_information) + _validate_hoomd_ff( + self.force_field.hoomd_forces, self.topology_information + ) elif isinstance(self.force_field, List): _validate_hoomd_ff(self.force_field, self.topology_information) else: diff --git a/flowermd/base/simulation.py b/flowermd/base/simulation.py index 2afcbc4d..04ba2d44 100644 --- a/flowermd/base/simulation.py +++ b/flowermd/base/simulation.py @@ -71,7 +71,8 @@ def __init__( ): if not isinstance(forcefield, Iterable) or isinstance(forcefield, str): raise ValueError( - "forcefield must be a sequence of " "hoomd.md.force.Force objects." + "forcefield must be a sequence of " + "hoomd.md.force.Force objects." ) else: for obj in forcefield: @@ -86,7 +87,9 @@ def __init__( self.gsd_write_freq = int(gsd_write_freq) self.maximum_write_buffer_size = gsd_max_buffer_size self.log_write_freq = int(log_write_freq) - self._std_out_freq = int((self.gsd_write_freq + self.log_write_freq) / 2) + self._std_out_freq = int( + (self.gsd_write_freq + self.log_write_freq) / 2 + ) self.gsd_file_name = gsd_file_name self.log_file_name = log_file_name self.log_quantities = [ @@ -433,7 +436,9 @@ def thermostat(self, thermostat): thermostat : flowermd.utils.HOOMDThermostats, required The type of thermostat to use. """ - if not issubclass(self._thermostat, hoomd.md.methods.thermostats.Thermostat): + if not issubclass( + self._thermostat, hoomd.md.methods.thermostats.Thermostat + ): raise ValueError( f"Invalid thermostat. Please choose from: {HOOMDThermostats}" ) @@ -522,7 +527,9 @@ def _initialize_thermostat(self, thermostat_kwargs): required_thermostat_kwargs = {} for k in inspect.signature(self.thermostat).parameters: if k not in thermostat_kwargs.keys(): - raise ValueError(f"Missing required parameter {k} for thermostat.") + raise ValueError( + f"Missing required parameter {k} for thermostat." + ) required_thermostat_kwargs[k] = thermostat_kwargs[k] return self.thermostat(**required_thermostat_kwargs) @@ -544,7 +551,9 @@ def set_integrator_method(self, integrator_method, method_kwargs): if not self.integrator: # Integrator and method not yet created self.integrator = hoomd.md.Integrator( dt=self.dt, - integrate_rotational_dof=(True if self._rigid_constraint else False), + integrate_rotational_dof=( + True if self._rigid_constraint else False + ), ) if self._rigid_constraint: self.integrator.rigid = self._rigid_constraint @@ -706,7 +715,9 @@ def run_update_volume( self.set_integrator_method( integrator_method=hoomd.md.methods.ConstantVolume, method_kwargs={ - "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), + "thermostat": self._initialize_thermostat( + {"kT": kT, "tau": tau_kt} + ), "filter": self.integrate_group, }, ) @@ -838,7 +849,9 @@ def run_NPT( "rescale_all": rescale_all, "gamma": gamma, "filter": self.integrate_group, - "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), + "thermostat": self._initialize_thermostat( + {"kT": kT, "tau": tau_kt} + ), }, ) if thermalize_particles: @@ -881,7 +894,9 @@ def run_NVT( self.set_integrator_method( integrator_method=hoomd.md.methods.ConstantVolume, method_kwargs={ - "thermostat": self._initialize_thermostat({"kT": kT, "tau": tau_kt}), + "thermostat": self._initialize_thermostat( + {"kT": kT, "tau": tau_kt} + ), "filter": self.integrate_group, }, ) @@ -1090,13 +1105,17 @@ def _thermalize_system(self, kT): filter=self.integrate_group, kT=kT.range[0] ) else: - self.state.thermalize_particle_momenta(filter=self.integrate_group, kT=kT) + self.state.thermalize_particle_momenta( + filter=self.integrate_group, kT=kT + ) def _lj_force(self): """Return the Lennard-Jones pair force.""" if not self.integrator: lj_force = [ - f for f in self._forcefield if isinstance(f, hoomd.md.pair.pair.LJ) + f + for f in self._forcefield + if isinstance(f, hoomd.md.pair.pair.LJ) ][0] else: lj_force = [ @@ -1122,7 +1141,9 @@ def _create_state(self, initial_state): print("Initializing simulation state from a GSD file.") self.create_state_from_gsd(initial_state) elif isinstance(initial_state, hoomd.snapshot.Snapshot): - print("Initializing simulation state from a hoomd.snapshot.Snapshot") + print( + "Initializing simulation state from a hoomd.snapshot.Snapshot" + ) self.create_state_from_snapshot(initial_state) elif isinstance(initial_state, gsd.hoomd.Frame): print("Initializing simulation state from a gsd.hoomd.Frame.") @@ -1130,7 +1151,9 @@ def _create_state(self, initial_state): def _add_hoomd_writers(self): """Create gsd and log writers.""" - gsd_logger = hoomd.logging.Logger(categories=["scalar", "string", "sequence"]) + gsd_logger = hoomd.logging.Logger( + categories=["scalar", "string", "sequence"] + ) logger = hoomd.logging.Logger(categories=["scalar", "string"]) gsd_logger.add(self, quantities=["timestep", "tps"]) logger.add(self, quantities=["timestep", "tps"]) diff --git a/flowermd/base/system.py b/flowermd/base/system.py index a9609160..6416be53 100644 --- a/flowermd/base/system.py +++ b/flowermd/base/system.py @@ -89,12 +89,16 @@ def __init__( if isinstance(mol_item, Molecule): # keep track of molecule types indices to assign to sites # before applying forcefield - self._mol_type_idx.extend([self.n_mol_types] * mol_item.n_particles) + self._mol_type_idx.extend( + [self.n_mol_types] * mol_item.n_particles + ) self.all_molecules.extend(mol_item.molecules) # if ff is provided in the Molecule class use that as the ff if mol_item.force_field: if isinstance(mol_item.force_field, BaseHOOMDForcefield): - self._hoomd_forcefield.extend(mol_item.force_field.hoomd_forces) + self._hoomd_forcefield.extend( + mol_item.force_field.hoomd_forces + ) elif isinstance(mol_item.force_field, BaseXMLForcefield): self._gmso_forcefields_dict[str(self.n_mol_types)] = ( mol_item.force_field.gmso_ff @@ -151,7 +155,9 @@ def mass(self): @property def net_charge(self): """Net charge of the system.""" - return sum(site.charge if site.charge else 0 for site in self.gmso_system.sites) + return sum( + site.charge if site.charge else 0 for site in self.gmso_system.sites + ) @property def box(self): @@ -303,8 +309,13 @@ def hoomd_snapshot(self): @property def hoomd_forcefield(self): """List of HOOMD forces.""" - if self._ff_refs != self.reference_values and self._gmso_forcefields_dict: - self._hoomd_forcefield = self._create_hoomd_forcefield(**self._ff_kwargs) + if ( + self._ff_refs != self.reference_values + and self._gmso_forcefields_dict + ): + self._hoomd_forcefield = self._create_hoomd_forcefield( + **self._ff_kwargs + ) return self._hoomd_forcefield def remove_hydrogens(self): @@ -316,7 +327,9 @@ def remove_hydrogens(self): """ # Try by element first: hydrogens = [ - site for site in self.gmso_system.sites if site.element.atomic_number == 1 + site + for site in self.gmso_system.sites + if site.element.atomic_number == 1 ] # If none found by element; try by mass if len(hydrogens) == 0: @@ -326,7 +339,9 @@ def remove_hydrogens(self): if site.mass.to("amu").value == 1.008 ] if len(hydrogens) == 0: - warnings.warn("Hydrogen atoms could not be found by element or mass") + warnings.warn( + "Hydrogen atoms could not be found by element or mass" + ) for h in hydrogens: # Find bond and other site in bond, add mass and charge for bond in self.gmso_system.iter_connections_by_site( @@ -350,15 +365,18 @@ def _scale_charges(self): """ charges = np.array( - [site.charge if site.charge else 0 for site in self.gmso_system.sites] + [ + site.charge if site.charge else 0 + for site in self.gmso_system.sites + ] ) net_charge = sum(charges) abs_charge = sum(abs(charges)) if abs_charge != 0: for site in self.gmso_system.sites: - site.charge -= abs(site.charge if site.charge else 0 * u.Unit("C")) * ( - net_charge / abs_charge - ) + site.charge -= abs( + site.charge if site.charge else 0 * u.Unit("C") + ) * (net_charge / abs_charge) def pickle_forcefield(self, file_path="forcefield.pickle"): """Pickle the list of HOOMD forces. @@ -401,7 +419,9 @@ def _create_hoomd_forcefield(self, r_cut, nlist_buffer, pppm_kwargs): nlist_buffer=nlist_buffer, pppm_kwargs=pppm_kwargs, auto_scale=False, - base_units=(self._reference_values if self._reference_values else None), + base_units=( + self._reference_values if self._reference_values else None + ), ) for force in ff: force_list.extend(ff[force]) @@ -413,7 +433,9 @@ def _create_hoomd_snapshot(self): snap, refs = to_gsd_snapshot( top=self.gmso_system, auto_scale=False, - base_units=(self._reference_values if self._reference_values else None), + base_units=( + self._reference_values if self._reference_values else None + ), ) self._snap_refs = self._reference_values.copy() return snap @@ -460,7 +482,9 @@ def _validate_forcefield(self, input_forcefield): else: # if there is only one ff for all molecule types ff_index = 0 - self._gmso_forcefields_dict[str(i)] = _force_field[ff_index].gmso_ff + self._gmso_forcefields_dict[str(i)] = _force_field[ + ff_index + ].gmso_ff def _assign_site_mol_type_idx(self): """Assign molecule type index to the gmso sites.""" @@ -539,17 +563,24 @@ def apply_forcefield( if not self._reference_values: epsilons = [ - s.atom_type.parameters["epsilon"] for s in self.gmso_system.sites + s.atom_type.parameters["epsilon"] + for s in self.gmso_system.sites + ] + sigmas = [ + s.atom_type.parameters["sigma"] for s in self.gmso_system.sites ] - sigmas = [s.atom_type.parameters["sigma"] for s in self.gmso_system.sites] masses = [s.mass for s in self.gmso_system.sites] energy_scale = np.max(epsilons) if self.auto_scale else 1.0 length_scale = np.max(sigmas) if self.auto_scale else 1.0 mass_scale = np.max(masses) if self.auto_scale else 1.0 - self._reference_values["energy"] = energy_scale * epsilons[0].unit_array - self._reference_values["length"] = length_scale * sigmas[0].unit_array + self._reference_values["energy"] = ( + energy_scale * epsilons[0].unit_array + ) + self._reference_values["length"] = ( + length_scale * sigmas[0].unit_array + ) self._reference_values["mass"] = mass_scale * masses[0].unit_array if remove_hydrogens: @@ -571,7 +602,9 @@ def visualize(self): if self.system: self.system.visualize().show() else: - raise ValueError("The initial configuraiton has not been created yet.") + raise ValueError( + "The initial configuraiton has not been created yet." + ) class Pack(System): @@ -627,7 +660,8 @@ def __init__( if not isinstance(density, u.array.unyt_quantity): self.density = density * u.Unit("g") / u.Unit("cm**3") warnings.warn( - "Units for density were not given, assuming " "units of g/cm**3." + "Units for density were not given, assuming " + "units of g/cm**3." ) else: self.density = density @@ -721,7 +755,9 @@ def _build_system(self): try: comp1 = self.all_molecules[next_idx] comp2 = self.all_molecules[next_idx + 1] - comp2.translate(self.basis_vector * np.array([self.x, self.y, 0])) + comp2.translate( + self.basis_vector * np.array([self.x, self.y, 0]) + ) # TODO: what if comp1 and comp2 have different names? unit_cell = mb.Compound( subcompounds=[comp1, comp2], name=comp1.name @@ -740,5 +776,7 @@ def _build_system(self): z_len = bounding_box.lengths[2] + 0.2 # Center the lattice in its box system.box = mb.box.Box(np.array([x_len, y_len, z_len])) - system.translate_to((system.box.Lx / 2, system.box.Ly / 2, system.box.Lz / 2)) + system.translate_to( + (system.box.Lx / 2, system.box.Ly / 2, system.box.Lz / 2) + ) return system diff --git a/flowermd/internal/ff_utils.py b/flowermd/internal/ff_utils.py index 5d5ab8b0..932a0054 100644 --- a/flowermd/internal/ff_utils.py +++ b/flowermd/internal/ff_utils.py @@ -20,7 +20,9 @@ def ff_xml_directory(): for dirpath, dirnames, filenames in os.walk(FF_DIR): for file in filenames: if file.endswith(".xml"): - ff_xml_directory[file.split(".xml")[0]] = os.path.join(dirpath, file) + ff_xml_directory[file.split(".xml")[0]] = os.path.join( + dirpath, file + ) return ff_xml_directory diff --git a/flowermd/internal/utils.py b/flowermd/internal/utils.py index 3ebb44b7..a8f6bea1 100644 --- a/flowermd/internal/utils.py +++ b/flowermd/internal/utils.py @@ -77,7 +77,9 @@ def _is_float(num): raise ValueError("The reference value is not a number.") # if ref_value is an instance of unyt_quantity, check the dimension. - if isinstance(ref_value, u.unyt_quantity) and _is_valid_dimension(ref_value.units): + if isinstance(ref_value, u.unyt_quantity) and _is_valid_dimension( + ref_value.units + ): return ref_value # if ref_value is a string, check if it is a number and if it is, check if # the unit exists in unyt and has the correct dimension. diff --git a/flowermd/library/forcefields.py b/flowermd/library/forcefields.py index 313218f7..3721ef8b 100644 --- a/flowermd/library/forcefields.py +++ b/flowermd/library/forcefields.py @@ -61,7 +61,9 @@ class OPLS_AA_DIMETHYLETHER(BaseXMLForcefield): """OPLS All Atom for dimethyl ether molecule forcefield class.""" def __init__(self, forcefield_files=f"{FF_DIR}/dimethylether_opls.xml"): - super(OPLS_AA_DIMETHYLETHER, self).__init__(forcefield_files=forcefield_files) + super(OPLS_AA_DIMETHYLETHER, self).__init__( + forcefield_files=forcefield_files + ) self.description = ( "Based on flowermd.forcefields.OPLS_AA. " "Trimmed down to include only dimethyl ether parameters." @@ -406,7 +408,9 @@ def _load_file(file, **kwargs): pair_dict[pair_type]["U"] = table[:, 1] pair_dict[pair_type]["F"] = table[:, 2] if len(pair_r_min) != len(pair_r_max) != 1: - raise ValueError("All pair files must have the same r-range values") + raise ValueError( + "All pair files must have the same r-range values" + ) # Read bond files bond_dict = dict() if bonds: @@ -426,9 +430,12 @@ def _load_file(file, **kwargs): for angle_type in angles: table = _load_file(angles[angle_type], **kwargs) thetas = table[:, 0] - if thetas[0] != 0 or not np.allclose(thetas[-1], np.pi, atol=1e-5): + if thetas[0] != 0 or not np.allclose( + thetas[-1], np.pi, atol=1e-5 + ): raise ValueError( - "Angle values must be evenly spaced and " "range from 0 to Pi." + "Angle values must be evenly spaced and " + "range from 0 to Pi." ) angle_dict[angle_type] = dict() angle_dict[angle_type]["U"] = table[:, 1] @@ -439,9 +446,9 @@ def _load_file(file, **kwargs): for dih_type in dihedrals: table = _load_file(dihedrals[dih_type], **kwargs) thetas = table[:, 0] - if not np.allclose(thetas[0], -np.pi, atol=1e-5) or not np.allclose( - thetas[-1], np.pi, atol=1e-5 - ): + if not np.allclose( + thetas[0], -np.pi, atol=1e-5 + ) or not np.allclose(thetas[-1], np.pi, atol=1e-5): raise ValueError( "Dihedral angle values must be evenly spaced and " "range from -Pi to Pi." @@ -467,7 +474,9 @@ def _create_forcefield(self): nlist = hoomd.md.nlist.Cell( buffer=self.nlist_buffer, exclusions=self.exclusions ) - pair_table = hoomd.md.pair.Table(nlist=nlist, default_r_cut=self.r_cut) + pair_table = hoomd.md.pair.Table( + nlist=nlist, default_r_cut=self.r_cut + ) for pair_type in self.pairs: U = self.pairs[pair_type]["U"] F = self.pairs[pair_type]["F"] @@ -475,7 +484,9 @@ def _create_forcefield(self): raise ValueError( "The energy and force arrays are not the same size." ) - pair_table.params[tuple(pair_type)] = dict(r_min=self.r_min, U=U, F=F) + pair_table.params[tuple(pair_type)] = dict( + r_min=self.r_min, U=U, F=F + ) forces.append(pair_table) # Create bond forces if self.bonds: diff --git a/flowermd/library/polymers.py b/flowermd/library/polymers.py index ac99d311..2ac9038f 100644 --- a/flowermd/library/polymers.py +++ b/flowermd/library/polymers.py @@ -258,7 +258,9 @@ def _build(self, length): bead_pair = "-".join([last_bead.name, next_bead.name]) bond_length = self.bond_lengths.get(bead_pair, None) if not bond_length: - bead_pair_rev = "-".join([next_bead.name, last_bead.name]) + bead_pair_rev = "-".join( + [next_bead.name, last_bead.name] + ) bond_length = self.bond_lengths.get(bead_pair_rev, None) if not bond_length: raise ValueError( @@ -317,8 +319,12 @@ def __init__(self, lengths, num_mols, lpar, bead_mass, bond_length): def _build(self, length): # Build up ellipsoid bead bead = mb.Compound(name="ellipsoid") - head = mb.Compound(pos=(0, 0, self.lpar), name="A", mass=self.bead_mass / 4) - tail = mb.Compound(pos=(0, 0, -self.lpar), name="A", mass=self.bead_mass / 4) + head = mb.Compound( + pos=(0, 0, self.lpar), name="A", mass=self.bead_mass / 4 + ) + tail = mb.Compound( + pos=(0, 0, -self.lpar), name="A", mass=self.bead_mass / 4 + ) head_mid = mb.Compound( pos=(0, 0, self.lpar / 2), name="B", mass=self.bead_mass / 4 ) diff --git a/flowermd/library/simulations/tensile.py b/flowermd/library/simulations/tensile.py index c74d5b3c..7b2b37ea 100644 --- a/flowermd/library/simulations/tensile.py +++ b/flowermd/library/simulations/tensile.py @@ -66,12 +66,16 @@ def __init__( self.fix_right = hoomd.filter.Tags(right_tags.astype(np.uint32)) all_fixed = hoomd.filter.Union(self.fix_left, self.fix_right) # Set the group of particles to be integrated over - self.integrate_group = hoomd.filter.SetDifference(hoomd.filter.All(), all_fixed) + self.integrate_group = hoomd.filter.SetDifference( + hoomd.filter.All(), all_fixed + ) @property def strain(self): """The current strain of the simulation.""" - delta_L = self.box_lengths_reduced[self._axis_index] - self.initial_length + delta_L = ( + self.box_lengths_reduced[self._axis_index] - self.initial_length + ) return delta_L / self.initial_length def run_tensile(self, strain, n_steps, kT, tau_kt, period): diff --git a/flowermd/modules/surface_wetting/surface_wetting.py b/flowermd/modules/surface_wetting/surface_wetting.py index e58f84fc..3f35383e 100644 --- a/flowermd/modules/surface_wetting/surface_wetting.py +++ b/flowermd/modules/surface_wetting/surface_wetting.py @@ -98,11 +98,12 @@ def run_droplet( """ # Shrink down to high density - if not isinstance(shrink_density, u.array.unyt_quantity) and not isinstance( - final_density, u.array.unyt_quantity - ): + if not isinstance( + shrink_density, u.array.unyt_quantity + ) and not isinstance(final_density, u.array.unyt_quantity): warnings.warn( - "Units for density were not given, assuming " "units of g/cm**3." + "Units for density were not given, assuming " + "units of g/cm**3." ) target_box_shrink = get_target_box_mass_density( density=shrink_density * (u.g / (u.cm**3)), @@ -122,7 +123,10 @@ def run_droplet( target_box_final = get_target_box_mass_density( density=final_density, mass=self.mass.to("g") ) - elif shrink_density.units.dimensions == number_density.units.dimensions: + elif ( + shrink_density.units.dimensions + == number_density.units.dimensions + ): raise ValueError( "For now, only mass density is supported " "in the surface wetting module." @@ -249,7 +253,9 @@ def _build_snapshot(self): ] self.drop_ptypes = self.drop_snapshot.particles.types - wetting_snapshot.particles.types = self.surface_ptypes + self.drop_ptypes + wetting_snapshot.particles.types = ( + self.surface_ptypes + self.drop_ptypes + ) wetting_snapshot.particles.typeid = np.concatenate( ( @@ -332,7 +338,8 @@ def _build_snapshot(self): self.surface_snapshot.dihedrals.N + self.drop_snapshot.dihedrals.N ) wetting_snapshot.dihedrals.types = ( - self.surface_snapshot.dihedrals.types + self.drop_snapshot.dihedrals.types + self.surface_snapshot.dihedrals.types + + self.drop_snapshot.dihedrals.types ) wetting_snapshot.dihedrals.typeid = np.concatenate( ( @@ -404,7 +411,9 @@ def _adjust_particle_positions(self, wetting_box): self.drop_snapshot.particles.position, axis=0 ) # shift drop particles z position to be at the top of surface - z_shift = np.abs(min(drop_pos[:, 2]) - max(surface_pos[:, 2])) - self.gap + z_shift = ( + np.abs(min(drop_pos[:, 2]) - max(surface_pos[:, 2])) - self.gap + ) drop_pos[:, 2] -= z_shift wetting_pos = np.concatenate((surface_pos, drop_pos), axis=0) return wetting_pos diff --git a/flowermd/modules/surface_wetting/utils.py b/flowermd/modules/surface_wetting/utils.py index 11a90dd9..c159530a 100644 --- a/flowermd/modules/surface_wetting/utils.py +++ b/flowermd/modules/surface_wetting/utils.py @@ -102,7 +102,9 @@ def _combine_lj_forces( surface_ptype, drop_ptype, ) not in list(lj.params.keys()): - drop_epsilon = drop_lj.params[(drop_ptype, drop_ptype)]["epsilon"] + drop_epsilon = drop_lj.params[(drop_ptype, drop_ptype)][ + "epsilon" + ] surface_epsilon = surface_lj.params[ (surface_ptype[1:], surface_ptype[1:]) ]["epsilon"] diff --git a/flowermd/modules/welding/utils.py b/flowermd/modules/welding/utils.py index 623f463e..029aa467 100644 --- a/flowermd/modules/welding/utils.py +++ b/flowermd/modules/welding/utils.py @@ -50,9 +50,13 @@ def add_void_particles( ) # Set updated mass and charges init_mass = snapshot.particles.mass - snapshot.particles.mass = np.concatenate((init_mass, np.array([1])), axis=None) + snapshot.particles.mass = np.concatenate( + (init_mass, np.array([1])), axis=None + ) init_charges = snapshot.particles.charge - snapshot.particles.charge = np.concatenate((init_charges, np.array([0])), axis=None) + snapshot.particles.charge = np.concatenate( + (init_charges, np.array([0])), axis=None + ) # Updated LJ params lj = [i for i in forcefield if isinstance(i, hoomd.md.pair.LJ)][0] for ptype in snapshot.particles.types: diff --git a/flowermd/modules/welding/welding.py b/flowermd/modules/welding/welding.py index f2518bc5..3894b63d 100644 --- a/flowermd/modules/welding/welding.py +++ b/flowermd/modules/welding/welding.py @@ -62,7 +62,9 @@ def _build(self): void_id = snap.particles.types.index("VOID") snap.particles.types.remove("VOID") keep_indices = np.where(snap.particles.typeid != void_id)[0] - snap.particles.position = snap.particles.position[keep_indices] + snap.particles.position = snap.particles.position[ + keep_indices + ] snap.particles.mass = snap.particles.mass[keep_indices] snap.particles.charge = snap.particles.charge[keep_indices] snap.particles.orientation = snap.particles.orientation[ @@ -87,14 +89,18 @@ def _build(self): # Set up snapshot.particles info: # Get set of new coordiantes, shifted along interface axis - shift = (snap_L.configuration.box[axis_index] + self.gap - self.wall_sigma) / 2 + shift = ( + snap_L.configuration.box[axis_index] + self.gap - self.wall_sigma + ) / 2 right_pos = np.copy(snap_R.particles.position) right_pos[:, axis_index] += shift left_pos = np.copy(snap_L.particles.position) left_pos[:, axis_index] -= shift pos = np.concatenate((left_pos, right_pos), axis=None) - mass = np.concatenate((snap_L.particles.mass, snap_R.particles.mass), axis=None) + mass = np.concatenate( + (snap_L.particles.mass, snap_R.particles.mass), axis=None + ) charges = np.concatenate( (snap_L.particles.charge, snap_R.particles.charge), axis=None ) @@ -110,7 +116,9 @@ def _build(self): # Set up bonds: bond_group_left = np.copy(snap_L.bonds.group) bond_group_right = np.copy(snap_R.bonds.group) + snap_R.particles.N - bond_group = np.concatenate((bond_group_left, bond_group_right), axis=None) + bond_group = np.concatenate( + (bond_group_left, bond_group_right), axis=None + ) bond_type_ids = np.concatenate( (snap_L.bonds.typeid, snap_R.bonds.typeid), axis=None ) @@ -121,7 +129,9 @@ def _build(self): # Set up angles: angle_group_left = np.copy(snap_L.angles.group) angle_group_right = np.copy(snap_R.angles.group) + snap_L.particles.N - angle_group = np.concatenate((angle_group_left, angle_group_right), axis=None) + angle_group = np.concatenate( + (angle_group_left, angle_group_right), axis=None + ) angle_type_ids = np.concatenate( (snap_L.angles.typeid, snap_R.angles.typeid), axis=None ) @@ -131,7 +141,9 @@ def _build(self): # Set up dihedrals: dihedral_group_left = np.copy(snap_L.dihedrals.group) - dihedral_group_right = np.copy(snap_R.dihedrals.group) + snap_L.particles.N + dihedral_group_right = ( + np.copy(snap_R.dihedrals.group) + snap_L.particles.N + ) dihedral_group = np.concatenate( (dihedral_group_left, dihedral_group_right), axis=None ) diff --git a/flowermd/tests/base/test_molecule.py b/flowermd/tests/base/test_molecule.py index d0ef0d5a..38c3b015 100644 --- a/flowermd/tests/base/test_molecule.py +++ b/flowermd/tests/base/test_molecule.py @@ -43,7 +43,9 @@ def test_molecule_topology_benzene(self, benzene_mb): assert not any(molecule.topology_information["particle_charge"]) def test_validate_force_field_oplsaa(self, benzene_mb): - molecule = Molecule(num_mols=2, force_field=OPLS_AA(), compound=benzene_mb) + molecule = Molecule( + num_mols=2, force_field=OPLS_AA(), compound=benzene_mb + ) assert molecule.gmso_molecule.is_typed() assert set(molecule.topology_information["particle_types"]) == { "opls_145", @@ -64,15 +66,23 @@ def test_validate_force_field_xml_file_path(self, benzene_mb, benzene_xml): } assert any(molecule.topology_information["particle_charge"]) - def test_validate_force_field_hoomd_ff_aa(self, benzene_mb, benzene_hoomd_ff): + def test_validate_force_field_hoomd_ff_aa( + self, benzene_mb, benzene_hoomd_ff + ): hoomd_ff = benzene_hoomd_ff(include_hydrogen=True) - molecule = Molecule(num_mols=2, force_field=hoomd_ff, compound=benzene_mb) + molecule = Molecule( + num_mols=2, force_field=hoomd_ff, compound=benzene_mb + ) assert molecule.force_field == hoomd_ff assert not molecule.gmso_molecule.is_typed() - def test_validate_fore_field_hoomd_ff_ua(self, benzene_mb, benzene_hoomd_ff): + def test_validate_fore_field_hoomd_ff_ua( + self, benzene_mb, benzene_hoomd_ff + ): hoomd_ff = benzene_hoomd_ff(include_hydrogen=False) - molecule = Molecule(num_mols=2, force_field=hoomd_ff, compound=benzene_mb) + molecule = Molecule( + num_mols=2, force_field=hoomd_ff, compound=benzene_mb + ) assert molecule.force_field == hoomd_ff assert not molecule.gmso_molecule.is_typed() @@ -155,7 +165,9 @@ def test_coarse_grain_with_multiple_single_beads( assert molecule.molecules[0].n_particles == 4 assert molecule.topology_information["bond_types"] == {("A", "A")} - def test_coarse_grain_with_different_beads(self, pps_smiles, benzene_smiles): + def test_coarse_grain_with_different_beads( + self, pps_smiles, benzene_smiles + ): molecule = Molecule(num_mols=2, smiles=pps_smiles) molecule.coarse_grain(beads={"A": benzene_smiles, "B": "S"}) assert molecule.molecules[0].n_particles == 2 @@ -181,7 +193,9 @@ def test_polymer(self, dimethylether_smiles): assert polymer.n_particles == 23 assert polymer.n_bonds == 22 assert ("O", "C", "C") in polymer.topology_information["angle_types"] - assert ("O", "C", "C", "O") in polymer.topology_information["dihedral_types"] + assert ("O", "C", "C", "O") in polymer.topology_information[ + "dihedral_types" + ] def test_polymer_different_chain_lengths(self, dimethylether_smiles): polymer = Polymer( @@ -233,7 +247,9 @@ def test_copolymer_with_sequence(self, polyethylene, polyDME): ) assert copolymer.n_particles == 22 assert copolymer.random_sequence is False - assert ("C", "C", "C", "C") in copolymer.topology_information["dihedral_types"] + assert ("C", "C", "C", "C") in copolymer.topology_information[ + "dihedral_types" + ] def test_copolymer_with_sequence_different_chain_lengths( self, polyethylene, polyDME @@ -249,7 +265,9 @@ def test_copolymer_with_sequence_different_chain_lengths( assert copolymer.molecules[0].n_particles == 42 assert copolymer.molecules[1].n_particles == 62 - def test_copolymer_with_sequence_different_num_mol(self, polyethylene, polyDME): + def test_copolymer_with_sequence_different_num_mol( + self, polyethylene, polyDME + ): copolymer = CoPolymer( monomer_A=polyDME, monomer_B=polyethylene, @@ -272,7 +290,9 @@ def test_copolymer_random_sequence(self, polyethylene, polyDME): ) # sequence is BAA assert copolymer.n_particles == 22 - assert ("O", "C", "C", "O") in copolymer.topology_information["dihedral_types"] + assert ("O", "C", "C", "O") in copolymer.topology_information[ + "dihedral_types" + ] def test_align_z_axis(self, polyethylene): pe = polyethylene(num_mols=5, lengths=20) diff --git a/flowermd/tests/base/test_simulation.py b/flowermd/tests/base/test_simulation.py index 3c4b5271..a756f808 100644 --- a/flowermd/tests/base/test_simulation.py +++ b/flowermd/tests/base/test_simulation.py @@ -27,7 +27,9 @@ def test_initialize_from_system(self, benzene_system): def test_initialize_from_system_separate_ff( self, benzene_cg_system, cg_single_bead_ff ): - sim = Simulation.from_system(benzene_cg_system, forcefield=cg_single_bead_ff) + sim = Simulation.from_system( + benzene_cg_system, forcefield=cg_single_bead_ff + ) sim.run_NVT(kT=0.1, tau_kt=10, n_steps=500) def test_initialize_from_system_missing_ff(self, benzene_cg_system): @@ -57,7 +59,9 @@ def test_reference_values(self, benzene_system): forcefield=benzene_system.hoomd_forcefield, reference_values=benzene_system.reference_values, ) - assert np.isclose(float(sim.mass.value), benzene_system.mass.value, atol=1e-4) + assert np.isclose( + float(sim.mass.value), benzene_system.mass.value, atol=1e-4 + ) assert np.allclose(benzene_system.box.lengths, sim.box_lengths.value) def test_set_ref_values(self, benzene_system): @@ -198,13 +202,17 @@ def test_update_volume_by_density_factor(self, benzene_system): period=1, final_box_lengths=target_box, ) - assert np.isclose(sim.density.value, (init_density * 5).value, atol=1e-4) + assert np.isclose( + sim.density.value, (init_density * 5).value, atol=1e-4 + ) def test_change_methods(self, benzene_system): sim = Simulation.from_system(benzene_system) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=0) assert isinstance(sim.method, hoomd.md.methods.ConstantVolume) - sim.run_NPT(kT=1.0, tau_kt=0.01, tau_pressure=0.1, pressure=0.001, n_steps=0) + sim.run_NPT( + kT=1.0, tau_kt=0.01, tau_pressure=0.1, pressure=0.001, n_steps=0 + ) assert isinstance(sim.method, hoomd.md.methods.ConstantPressure) def test_change_dt(self, benzene_system): @@ -396,7 +404,9 @@ def test_gsd_logger(self, benzene_system): def test_bad_ff(self, benzene_system): with pytest.raises(ValueError): - Simulation(initial_state=benzene_system.hoomd_snapshot, forcefield="gaff") + Simulation( + initial_state=benzene_system.hoomd_snapshot, forcefield="gaff" + ) with pytest.raises(ValueError): Simulation( initial_state=benzene_system.hoomd_snapshot, diff --git a/flowermd/tests/base/test_system.py b/flowermd/tests/base/test_system.py index 5bb38de8..cade5f96 100644 --- a/flowermd/tests/base/test_system.py +++ b/flowermd/tests/base/test_system.py @@ -23,7 +23,9 @@ class TestSystem(BaseTest): def test_single_mol_type(self, benzene_molecule): benzene_mols = benzene_molecule(n_mols=3) system = Pack(molecules=[benzene_mols], density=0.8) - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(benzene_mols.molecules) assert system.gmso_system.is_typed() @@ -36,7 +38,9 @@ def test_multiple_mol_types(self, benzene_molecule, ethane_molecule): benzene_mol = benzene_molecule(n_mols=3) ethane_mol = ethane_molecule(n_mols=2) system = Pack(molecules=[benzene_mol, ethane_mol], density=0.8) - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) assert system.n_mol_types == 2 assert len(system.all_molecules) == len(benzene_mol.molecules) + len( ethane_mol.molecules @@ -65,9 +69,9 @@ def test_multiple_mol_types_different_ff( auto_scale=True, ) assert system.n_mol_types == 2 - assert len(system.all_molecules) == len(dimethylether_mol.molecules) + len( - pps_mol.molecules - ) + assert len(system.all_molecules) == len( + dimethylether_mol.molecules + ) + len(pps_mol.molecules) assert system.gmso_system.sites[0].group == "0" assert system.gmso_system.sites[-1].group == "1" assert system.gmso_system.is_typed() @@ -90,7 +94,9 @@ def test_multiple_mol_types_different_ff( def test_system_from_mol2_mol_parameterization(self, benzene_molecule_mol2): benzene_mol = benzene_molecule_mol2(n_mols=3) system = Pack(molecules=[benzene_mol], density=0.8) - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) assert system.gmso_system.is_typed() assert len(system.hoomd_forcefield) > 0 assert system.n_particles == system.hoomd_snapshot.particles.N @@ -104,7 +110,9 @@ def test_remove_hydrogen(self, benzene_molecule): remove_hydrogens=True, auto_scale=True, ) - assert not any([s.element.atomic_number == 1 for s in system.gmso_system.sites]) + assert not any( + [s.element.atomic_number == 1 for s in system.gmso_system.sites] + ) assert system.gmso_system.is_typed() assert len(system.hoomd_forcefield) > 0 assert list(system.hoomd_forcefield[0].params.keys()) == [ @@ -155,7 +163,9 @@ def test_remove_hydrogen_no_hydrogen(self, benzene_molecule): auto_scale=True, ) hydrogens = [ - site for site in system.gmso_system.sites if site.element.atomic_number == 1 + site + for site in system.gmso_system.sites + if site.element.atomic_number == 1 ] for h_site in hydrogens: system.gmso_system.remove_site(h_site) @@ -178,7 +188,9 @@ def test_add_mass_charges(self, benzene_molecule): assert site.charge == 0 snap = system.hoomd_snapshot - assert np.allclose(sum(snap.particles.mass), 6 * (12.011 + 1.008), atol=1e-4) + assert np.allclose( + sum(snap.particles.mass), 6 * (12.011 + 1.008), atol=1e-4 + ) assert sum(snap.particles.charge) == 0 def test_pack_box(self, benzene_molecule): @@ -201,9 +213,13 @@ def test_mass(self, pps_molecule): def test_ref_length(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) - assert np.allclose(system.reference_length.to("angstrom").value, 3.5, atol=1e-3) + assert np.allclose( + system.reference_length.to("angstrom").value, 3.5, atol=1e-3 + ) reduced_box = system.hoomd_snapshot.configuration.box[0:3] calc_box = reduced_box * system.reference_length.to("nm").value assert np.allclose(calc_box[0], system.box.Lx, atol=1e-2) @@ -213,7 +229,9 @@ def test_ref_length(self, polyethylene): def test_ref_mass(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) total_red_mass = sum(system.hoomd_snapshot.particles.mass) assert np.allclose( system.mass.value, @@ -224,7 +242,9 @@ def test_ref_mass(self, polyethylene): def test_ref_energy(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=5) system = Pack(molecules=[polyethylene], density=1.0, base_units=dict()) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) assert np.allclose( system.reference_energy.to("kcal/mol").value, 0.066, atol=1e-3 ) @@ -235,7 +255,9 @@ def test_rebuild_snapshot(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) assert system._snap_refs == system.reference_values assert system._ff_refs == system.reference_values init_snap = system.hoomd_snapshot @@ -255,7 +277,9 @@ def test_ref_values_no_autoscale(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_length = 1 * u.angstrom system.reference_energy = 1 * u.kcal / u.mol system.reference_mass = 1 * u.amu @@ -265,12 +289,16 @@ def test_ref_values_no_autoscale(self, polyethylene): ) assert dict(system.hoomd_forcefield[3].params)["opls_135", "opls_135"][ "epsilon" - ] == system.gmso_system.sites[0].atom_type.parameters["epsilon"].to("kcal/mol") + ] == system.gmso_system.sites[0].atom_type.parameters["epsilon"].to( + "kcal/mol" + ) def test_set_ref_values(self, polyethylene): polyethylene = polyethylene(lengths=5, num_mols=1) system = Pack(molecules=[polyethylene], density=1.0) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -287,7 +315,9 @@ def test_set_ref_values_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) ref_value_dict = { "length": "1 angstrom", "energy": "3 kcal/mol", @@ -304,7 +334,9 @@ def test_set_ref_values_missing_key(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -318,7 +350,9 @@ def test_set_ref_values_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -333,7 +367,9 @@ def test_set_ref_values_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) ref_value_dict = { "length": 1 * u.angstrom, "energy": 3.0 * u.kcal / u.mol, @@ -348,7 +384,9 @@ def test_set_ref_length(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_length = 1 * u.angstrom assert system.reference_length == 1 * u.angstrom @@ -358,7 +396,9 @@ def test_set_ref_length_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_length = 1.0 @@ -368,7 +408,9 @@ def test_ref_length_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_length = "1 angstrom" assert system.reference_length == 1 * u.angstrom @@ -378,7 +420,9 @@ def test_ref_length_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0" @@ -388,7 +432,9 @@ def test_ref_length_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0 invalid_unit" @@ -398,7 +444,9 @@ def test_ref_length_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_length = 1.0 * u.g @@ -408,7 +456,9 @@ def test_ref_length_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_length = "1.0 g" @@ -418,7 +468,9 @@ def test_ref_length_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) system.reference_length = 1 * u.angstrom assert system.reference_length == 1 * u.angstrom @@ -428,7 +480,9 @@ def test_set_ref_energy(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_energy = 1 * u.kcal / u.mol assert system.reference_energy == 1 * u.kcal / u.mol @@ -438,7 +492,9 @@ def test_set_ref_energy_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 @@ -448,7 +504,9 @@ def test_ref_energy_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_energy = "1 kJ" assert system.reference_energy == 1 * u.kJ @@ -458,7 +516,9 @@ def test_ref_energy_string_comb_units(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_energy = "1 kcal/mol" assert system.reference_energy == 1 * u.kcal / u.mol @@ -468,7 +528,9 @@ def test_ref_energy_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0" @@ -478,7 +540,9 @@ def test_ref_energy_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0 invalid_unit" @@ -488,7 +552,9 @@ def test_ref_energy_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 * u.g @@ -498,7 +564,9 @@ def test_ref_energy_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = "1.0 m" @@ -508,7 +576,9 @@ def test_set_ref_energy_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) system.reference_energy = 1 * u.kcal / u.mol assert system.reference_energy == 1 * u.kcal / u.mol @@ -518,7 +588,9 @@ def test_set_ref_mass(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_mass = 1.0 * u.amu assert system.reference_mass == 1.0 * u.amu @@ -529,7 +601,9 @@ def test_set_ref_mass_invalid_type(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_mass = 1.0 @@ -539,7 +613,9 @@ def test_ref_mass_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.reference_mass = "1 g" assert system.reference_mass == 1.0 * u.g @@ -549,7 +625,9 @@ def test_ref_mass_invalid_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0" @@ -559,7 +637,9 @@ def test_ref_mass_invalid_unit_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0 invalid_unit" @@ -569,7 +649,9 @@ def test_ref_mass_invalid_dimension(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_energy = 1.0 * u.m @@ -579,7 +661,9 @@ def test_ref_mass_invalid_dimension_string(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) with pytest.raises(ReferenceUnitError): system.reference_mass = "1.0 m" @@ -589,7 +673,9 @@ def test_set_ref_mass_auto_scale_true(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) system.reference_mass = 1.0 * u.amu assert system.reference_mass == 1.0 * u.amu @@ -601,7 +687,9 @@ def test_apply_forcefield_no_forcefield(self, polyethylene): density=1.0, ) with pytest.raises(ForceFieldError): - system.apply_forcefield(r_cut=2.5, force_field=None, auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=None, auto_scale=False + ) def test_apply_forcefield_no_forcefield_w_mol_ff(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=3, force_field=OPLS_AA()) @@ -621,7 +709,9 @@ def test_apply_forcefield_w_mol_ff(self, benzene_molecule): density=1.0, ) with pytest.raises(ForceFieldError): - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) def test_validate_forcefield_invalid_ff_type(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=1) @@ -672,7 +762,9 @@ def test_forcefield_list_hoomd_ff(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) hoomd_ff = system.hoomd_forcefield @@ -691,7 +783,9 @@ def test_lattice_polymer(self, polyethylene): y=1, n=4, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=True + ) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(polyethylene.molecules) @@ -701,7 +795,9 @@ def test_lattice_polymer(self, polyethylene): for mol_class in system._molecules: for mol in mol_class.molecules: backbone = get_backbone_vector(mol.xyz) - assert np.allclose(np.abs(backbone), np.array([0, 0, 1]), atol=1e-1) + assert np.allclose( + np.abs(backbone), np.array([0, 0, 1]), atol=1e-1 + ) def test_lattice_molecule(self, benzene_molecule): benzene_mol = benzene_molecule(n_mols=32) @@ -711,7 +807,9 @@ def test_lattice_molecule(self, benzene_molecule): y=1, n=4, ) - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) assert system.n_mol_types == 1 assert len(system.all_molecules) == len(benzene_mol.molecules) assert len(system.hoomd_forcefield) > 0 @@ -754,7 +852,9 @@ def test_to_gsd(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.to_gsd("test.gsd") assert os.path.isfile(os.path.join(os.getcwd(), "test.gsd")) os.remove(os.path.join(os.getcwd(), "test.gsd")) @@ -775,7 +875,9 @@ def test_pickle_ff(self, polyethylene): molecules=[polyethylene], density=1.0, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA()], auto_scale=False + ) system.pickle_forcefield("forcefield.pickle") assert os.path.isfile(os.path.join(os.getcwd(), "forcefield.pickle")) os.remove(os.path.join(os.getcwd(), "forcefield.pickle")) diff --git a/flowermd/tests/base_test.py b/flowermd/tests/base_test.py index 3b8b7867..3a70949d 100644 --- a/flowermd/tests/base_test.py +++ b/flowermd/tests/base_test.py @@ -103,7 +103,9 @@ def _hoomd_ff(include_hydrogen, invalid_pair=False): ) bonds = self.benzene_hoomd_bond(include_hydrogen=include_hydrogen) angles = self.benzene_hoomd_angle(include_hydrogen=include_hydrogen) - dihedrals = self.benzene_hoomd_dihedral(include_hydrogen=include_hydrogen) + dihedrals = self.benzene_hoomd_dihedral( + include_hydrogen=include_hydrogen + ) return [pairs, bonds, angles, dihedrals] return _hoomd_ff @@ -145,7 +147,9 @@ def _pps_molecule(n_mols): @pytest.fixture() def dimethylether_molecule(self, dimethylether_smiles): def _dimethylether_molecule(n_mols): - dimethylether = Molecule(num_mols=n_mols, smiles=dimethylether_smiles) + dimethylether = Molecule( + num_mols=n_mols, smiles=dimethylether_smiles + ) return dimethylether return _dimethylether_molecule @@ -217,7 +221,9 @@ def benzene_system(self, benzene_mb): molecules=[benzene], density=0.2 * u.g / u.cm**3, ) - system.apply_forcefield(r_cut=2.5, force_field=OPLS_AA(), auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=OPLS_AA(), auto_scale=True + ) return system @pytest.fixture() diff --git a/flowermd/tests/library/test_forcefields.py b/flowermd/tests/library/test_forcefields.py index fb927109..30ae0920 100644 --- a/flowermd/tests/library/test_forcefields.py +++ b/flowermd/tests/library/test_forcefields.py @@ -117,9 +117,15 @@ def test_ellipsoid_ff(self): angle_theta0=1.57, ) assert len(ellipsoid_ff.hoomd_forces) == 3 - assert isinstance(ellipsoid_ff.hoomd_forces[-1], hoomd.md.pair.aniso.GayBerne) - assert ("R", "R") in list(dict(ellipsoid_ff.hoomd_forces[-1].params).keys()) - assert ("B", "R") in list(dict(ellipsoid_ff.hoomd_forces[-1].params).keys()) + assert isinstance( + ellipsoid_ff.hoomd_forces[-1], hoomd.md.pair.aniso.GayBerne + ) + assert ("R", "R") in list( + dict(ellipsoid_ff.hoomd_forces[-1].params).keys() + ) + assert ("B", "R") in list( + dict(ellipsoid_ff.hoomd_forces[-1].params).keys() + ) assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["epsilon"] == 0.0 assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["lperp"] == 0.0 assert ellipsoid_ff.hoomd_forces[-1].params["A", "A"]["lpar"] == 0.0 @@ -149,7 +155,9 @@ def test_from_txt_file(self): def test_from_csv_file(self): pair_file = os.path.join(ASSETS_DIR, "lj_pair_table.csv") - ff = TableForcefield.from_files(pairs={("A", "A"): pair_file}, delimiter=",") + ff = TableForcefield.from_files( + pairs={("A", "A"): pair_file}, delimiter="," + ) pair_data = np.genfromtxt(pair_file, delimiter=",") assert ff.r_min == pair_data[:, 0][0] assert ff.r_cut == pair_data[:, 0][-1] diff --git a/flowermd/tests/library/test_tensile.py b/flowermd/tests/library/test_tensile.py index 7cedc9ec..161dfc38 100644 --- a/flowermd/tests/library/test_tensile.py +++ b/flowermd/tests/library/test_tensile.py @@ -14,7 +14,9 @@ def test_tensile(self): y=1.2, n=4, ) - system.apply_forcefield(r_cut=2.5, force_field=[OPLS_AA_PPS()], auto_scale=True) + system.apply_forcefield( + r_cut=2.5, force_field=[OPLS_AA_PPS()], auto_scale=True + ) tensile_sim = Tensile( initial_state=system.hoomd_snapshot, diff --git a/flowermd/tests/modules/test_surface_wetting.py b/flowermd/tests/modules/test_surface_wetting.py index 345beda9..b2dc908b 100644 --- a/flowermd/tests/modules/test_surface_wetting.py +++ b/flowermd/tests/modules/test_surface_wetting.py @@ -43,7 +43,9 @@ def test_droplet_sim(self, polyethylene): final_density=0.05 * u.g / u.cm**3, tau_kt=drop_sim.dt * 100, ) - assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) + assert np.isclose( + drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 + ) def test_droplet_sim_no_units(self, polyethylene): drop_mol = polyethylene(num_mols=10, lengths=5) @@ -74,7 +76,9 @@ def test_droplet_sim_no_units(self, polyethylene): final_density=0.05, tau_kt=drop_sim.dt * 100, ) - assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) + assert np.isclose( + drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 + ) def test_droplet_sim_bad_units(self, polyethylene): drop_mol = polyethylene(num_mols=10, lengths=5) @@ -105,7 +109,9 @@ def test_droplet_sim_bad_units(self, polyethylene): final_density=0.05 * (u.cm**-3), tau_kt=drop_sim.dt * 100, ) - assert np.isclose(drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2) + assert np.isclose( + drop_sim.density.to(u.g / u.cm**3).value, 0.05, atol=1e-2 + ) class TestInterfaceBuilder(BaseTest): @@ -154,7 +160,8 @@ def test_interface_builder( + drop_snapshot.particles.types ) assert np.isclose( - interface.hoomd_snapshot.configuration.box[2] * drop_refs["length"].value, + interface.hoomd_snapshot.configuration.box[2] + * drop_refs["length"].value, 15, atol=1e-2, ) @@ -182,7 +189,9 @@ def test_interface_builder( class TestWettingSimulation(BaseTest): - def test_wetting_sim(self, surface_wetting_init_snapshot, surface_wetting_init_ff): + def test_wetting_sim( + self, surface_wetting_init_snapshot, surface_wetting_init_ff + ): # load surface wetting init snapshot snapshot = gsd.hoomd.open(surface_wetting_init_snapshot)[0] # load ff from pickle diff --git a/flowermd/tests/modules/test_weld.py b/flowermd/tests/modules/test_weld.py index 06f87ec7..09ed0e42 100644 --- a/flowermd/tests/modules/test_weld.py +++ b/flowermd/tests/modules/test_weld.py @@ -90,7 +90,9 @@ def test_slab_sim_xaxis(self, polyethylene_system): forcefield=polyethylene_system.hoomd_forcefield, log_write_freq=2000, ) - assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) + assert any( + [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] + ) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_slab_sim_yaxis(self, polyethylene_system): @@ -100,7 +102,9 @@ def test_slab_sim_yaxis(self, polyethylene_system): interface_axis=(0, 1, 0), log_write_freq=2000, ) - assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) + assert any( + [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] + ) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_slab_sim_zaxis(self, polyethylene_system): @@ -110,7 +114,9 @@ def test_slab_sim_zaxis(self, polyethylene_system): interface_axis=(0, 0, 1), log_write_freq=2000, ) - assert any([isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces]) + assert any( + [isinstance(i, hoomd.md.external.wall.LJ) for i in sim.forces] + ) sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) def test_weld_sim(self, polyethylene_system): @@ -122,7 +128,9 @@ def test_weld_sim(self, polyethylene_system): sim.run_NVT(kT=1.0, tau_kt=0.01, n_steps=500) sim.save_restart_gsd() # Create interface system from slab restart.gsd file - interface = Interface(gsd_files="restart.gsd", interface_axis="x", gap=0.1) + interface = Interface( + gsd_files="restart.gsd", interface_axis="x", gap=0.1 + ) sim = WeldSimulation( initial_state=interface.hoomd_snapshot, forcefield=polyethylene_system.hoomd_forcefield, diff --git a/flowermd/tests/utils/test_actions.py b/flowermd/tests/utils/test_actions.py index 7cb2f90d..3955bfc3 100644 --- a/flowermd/tests/utils/test_actions.py +++ b/flowermd/tests/utils/test_actions.py @@ -11,18 +11,24 @@ class TestActions(BaseTest): def test_scale_epsilon(self, benzene_simulation): sim = benzene_simulation epsilon_scale = ScaleEpsilon(sim=sim, scale_factor=0.5) - energy_operation = hoomd.update.CustomUpdater(action=epsilon_scale, trigger=10) + energy_operation = hoomd.update.CustomUpdater( + action=epsilon_scale, trigger=10 + ) sim.operations.updaters.append(energy_operation) old_lj_force = copy.deepcopy(sim._lj_force().params) sim.run_NVT(n_steps=10, kT=1.0, tau_kt=1.0) new_lj_force = sim._lj_force().params for k in old_lj_force.keys(): - assert new_lj_force[k]["epsilon"] == old_lj_force[k]["epsilon"] + 0.5 + assert ( + new_lj_force[k]["epsilon"] == old_lj_force[k]["epsilon"] + 0.5 + ) def test_scale_sigma(self, benzene_simulation): sim = benzene_simulation sigma_scale = ScaleSigma(sim=sim, scale_factor=0.5) - energy_operation = hoomd.update.CustomUpdater(action=sigma_scale, trigger=10) + energy_operation = hoomd.update.CustomUpdater( + action=sigma_scale, trigger=10 + ) sim.operations.updaters.append(energy_operation) old_lj_force = copy.deepcopy(sim._lj_force().params) sim.run_NVT(n_steps=10, kT=1.0, tau_kt=1.0) diff --git a/flowermd/tests/utils/test_rigid_body.py b/flowermd/tests/utils/test_rigid_body.py index 9c5f5c52..08ed20d5 100644 --- a/flowermd/tests/utils/test_rigid_body.py +++ b/flowermd/tests/utils/test_rigid_body.py @@ -29,7 +29,8 @@ def test_ellipsoid_create_rigid_body(self): ) assert rigid_frame.particles.N == 8 + ellipsoid_chain.n_particles assert ( - rigid_frame.particles.types == ["R"] + system.hoomd_snapshot.particles.types + rigid_frame.particles.types + == ["R"] + system.hoomd_snapshot.particles.types ) assert rigid_frame.particles.mass[0] == 100 assert np.all( diff --git a/flowermd/tests/utils/test_utils.py b/flowermd/tests/utils/test_utils.py index 06283edd..505eb0d3 100644 --- a/flowermd/tests/utils/test_utils.py +++ b/flowermd/tests/utils/test_utils.py @@ -30,7 +30,9 @@ def test_validate_ref_value(self): assert validate_ref_value("1.0 kcal/mol", u.dimensions.energy) == ( 1.0 * u.kcal / u.mol ) - assert validate_ref_value("1.0 amu", u.dimensions.mass) == 1.0 * u.Unit("amu") + assert validate_ref_value("1.0 amu", u.dimensions.mass) == 1.0 * u.Unit( + "amu" + ) with pytest.raises(ReferenceUnitError): validate_ref_value("1.0 g", u.dimensions.energy) @@ -79,7 +81,9 @@ def test_target_box_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - target_box = get_target_box_number_density(density=density, n_beads=n_beads) + target_box = get_target_box_number_density( + density=density, n_beads=n_beads + ) L = target_box[0].value assert np.allclose(L**3, 100, atol=1e-8) @@ -87,7 +91,9 @@ def test_target_box_one_constraint_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - cubic_box = get_target_box_number_density(density=density, n_beads=n_beads) + cubic_box = get_target_box_number_density( + density=density, n_beads=n_beads + ) tetragonal_box = get_target_box_number_density( density=density, n_beads=n_beads, @@ -100,7 +106,9 @@ def test_target_box_two_constraint_number_density(self): sigma = 1 * u.nm n_beads = 100 density = 1 / sigma**3 - cubic_box = get_target_box_number_density(density=density, n_beads=n_beads) + cubic_box = get_target_box_number_density( + density=density, n_beads=n_beads + ) ortho_box = get_target_box_number_density( density=density, n_beads=n_beads, @@ -108,7 +116,9 @@ def test_target_box_two_constraint_number_density(self): y_constraint=cubic_box[0] / 2, ) assert cubic_box[0] / 2 == ortho_box[0] == ortho_box[1] != ortho_box[2] - assert np.allclose(ortho_box[2].value, cubic_box[0].value * 4, atol=1e-5) + assert np.allclose( + ortho_box[2].value, cubic_box[0].value * 4, atol=1e-5 + ) def test_calculate_box_length_bad_args(self): mass_density = 1 * u.g / (u.cm**3) @@ -122,12 +132,16 @@ def test_calculate_box_length_fixed_l_1d(self): mass = u.unyt_quantity(6.0, u.g) density = u.unyt_quantity(0.5, u.g / u.cm**3) fixed_L = u.unyt_quantity(3.0, u.cm) - box_length = _calculate_box_length(mass=mass, density=density, fixed_L=fixed_L) + box_length = _calculate_box_length( + mass=mass, density=density, fixed_L=fixed_L + ) assert box_length == 2.0 * u.cm def test_calculate_box_length_fixed_l_2d(self): mass = u.unyt_quantity(12.0, u.g) density = u.unyt_quantity(0.5, u.g / u.cm**3) fixed_L = u.unyt_array([3.0, 2.0], u.cm) - box_length = _calculate_box_length(mass=mass, density=density, fixed_L=fixed_L) + box_length = _calculate_box_length( + mass=mass, density=density, fixed_L=fixed_L + ) assert box_length == 4.0 * u.cm diff --git a/flowermd/utils/rigid_body.py b/flowermd/utils/rigid_body.py index fba8bb36..4484b1e6 100644 --- a/flowermd/utils/rigid_body.py +++ b/flowermd/utils/rigid_body.py @@ -70,7 +70,9 @@ def create_rigid_body( """ # find typeid sequence of the constituent particles types in a rigid bead p_types = np.array(snapshot.particles.types) - constituent_type_ids = np.where(p_types[:, None] == bead_constituents_types)[0] + constituent_type_ids = np.where( + p_types[:, None] == bead_constituents_types + )[0] # find indices that matches the constituent particle types typeids = snapshot.particles.typeid @@ -86,7 +88,9 @@ def create_rigid_body( n_rigid = rigid_const_idx.shape[0] # number of rigid bodies # calculate center of mass and its position for each rigid body - com_mass, com_position, com_moi = _get_com_mass_pos_moi(snapshot, rigid_const_idx) + com_mass, com_position, com_moi = _get_com_mass_pos_moi( + snapshot, rigid_const_idx + ) rigid_frame = gsd.hoomd.Frame() rigid_frame.particles.types = [bead_name] + snapshot.particles.types @@ -94,7 +98,9 @@ def create_rigid_body( rigid_frame.particles.typeid = np.concatenate( (([0] * n_rigid), snapshot.particles.typeid + 1) ) - rigid_frame.particles.mass = np.concatenate((com_mass, snapshot.particles.mass)) + rigid_frame.particles.mass = np.concatenate( + (com_mass, snapshot.particles.mass) + ) rigid_frame.particles.position = np.concatenate( (com_position, snapshot.particles.position) ) @@ -140,7 +146,9 @@ def create_rigid_body( # find local coordinates of the particles in the first rigid body # only need to find the local coordinates for the first rigid body - local_coords = snapshot.particles.position[rigid_const_idx[0]] - com_position[0] + local_coords = ( + snapshot.particles.position[rigid_const_idx[0]] - com_position[0] + ) rigid_constrain = hoomd.md.constrain.Rigid() rigid_constrain.body["R"] = { diff --git a/flowermd/utils/utils.py b/flowermd/utils/utils.py index e2c973ab..afb08949 100644 --- a/flowermd/utils/utils.py +++ b/flowermd/utils/utils.py @@ -88,7 +88,9 @@ def get_target_box_number_density( else: constraints = np.array([x_constraint, y_constraint, z_constraint]) fixed_L = constraints[np.not_equal(constraints, None).nonzero()] - L = _calculate_box_length(density=density, n_beads=n_beads, fixed_L=fixed_L) + L = _calculate_box_length( + density=density, n_beads=n_beads, fixed_L=fixed_L + ) constraints[np.equal(constraints, None).nonzero()] = L Lx, Ly, Lz = constraints return np.array([Lx, Ly, Lz]) * u.Unit(Lx.units) diff --git a/setup.py b/setup.py index 1b717a86..cd3f223a 100644 --- a/setup.py +++ b/setup.py @@ -6,7 +6,9 @@ # Package meta-data. NAME = "flowermd" -DESCRIPTION = "Package making it easier to build and simulate polymers in Hoomd-Blue." +DESCRIPTION = ( + "Package making it easier to build and simulate polymers in Hoomd-Blue." +) URL = "https://github.com/cmelab/flowerMD" EMAIL = "chrisjones4@u.boisestate.edu" AUTHOR = "CME Lab" @@ -56,7 +58,9 @@ def run(self): pass self.status("Building Source and Wheel (universal) distribution…") - os.system("{0} setup.py sdist bdist_wheel --universal".format(sys.executable)) + os.system( + "{0} setup.py sdist bdist_wheel --universal".format(sys.executable) + ) self.status("Uploading the package to PyPi via Twine…") os.system("twine upload dist/*")