diff --git a/CHANGELOG.md b/CHANGELOG.md index 430332c5..14b16d64 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -15,6 +15,8 @@ This changelog is effective from the 2025 releases. * MultiJob now supports generic Job types for the self.children attribute ### Changed +* `Molecule.readmol2` can now read non-integer bond orders +* `Molecule.writepdb` now writes default values for res, resnum, fix, and occ columns * `view` function uses stdin mode for AMSview, reducing overhead for image creation * `vasp_output_to_ams` now reads the MD time step from `POTIM` in the OUTCAR and labels molecular-dynamics runs (`IBRION = 0`) as such, instead of using a fixed default time step diff --git a/src/scm/plams/mol/molecule.py b/src/scm/plams/mol/molecule.py index d89fe6eb..aa22bfa2 100644 --- a/src/scm/plams/mol/molecule.py +++ b/src/scm/plams/mol/molecule.py @@ -3184,6 +3184,8 @@ def readmol2(self, f: IO, **other: Any) -> None: spl = line.split() if len(spl) < 4: raise FileError(f"readmol2: Error in {f.name} line {i+1}: not enough values in line") + if "." in spl[3]: + spl[3] = f"{int(float(spl[3]))}" try: atom1 = self.atoms[int(spl[1]) - 1] atom2 = self.atoms[int(spl[2]) - 1] @@ -3278,8 +3280,14 @@ def writepdb(self, f: IO, **other: Any) -> None: pdb = PDBHandler() for i, at in enumerate(self.atoms): pdbatom = PDBAtom() + pdbatom.name = f"{at.symbol.upper():<2}" + pdbatom.name = "%-2s" % (at.symbol.upper()) pdbatom.coords = at.coords pdbatom.element = at.symbol.upper() + pdbatom.res = "LIG" + pdbatom.resnum = "0" + pdbatom.occ = "1.00" + pdbatom.fix = "0.00" if "pdb" in at.properties: if "res" in at.properties.pdb: pdbatom.res = at.properties.pdb.res @@ -3287,6 +3295,10 @@ def writepdb(self, f: IO, **other: Any) -> None: pdbatom.resnum = at.properties.pdb.resnum if "name" in at.properties.pdb: pdbatom.name = at.properties.pdb.name + if "occ" in at.properties.pdb: + pdbatom.occ = at.properties.pdb.occ + if "fix" in at.properties.pdb: + pdbatom.fix = at.properties.pdb.fix pdb.add_atom(pdbatom) if len(self.lattice) > 0: pdb.set_lattice(self.lattice) diff --git a/unit_tests/test_molecule.py b/unit_tests/test_molecule.py index 66a3e402..41479000 100644 --- a/unit_tests/test_molecule.py +++ b/unit_tests/test_molecule.py @@ -1,4 +1,5 @@ import os +from io import StringIO import pytest from abc import ABC, abstractmethod @@ -1113,6 +1114,42 @@ def test_write_multiple_molecules_to_pdb(pdb_folder, tmp_path): np.testing.assert_allclose(at.coords, at_ref.coords, atol=1e-08) +def test_simple_pdb_writing(xyz_folder, tmp_path): + """ + Test that a reasonable PDB format is written for a simple benzene molecule + + The file should contain for each atom a default resname and resnum, + as well as default fix and occ columns. + """ + benzene = Molecule(os.path.join(xyz_folder, "benzene.xyz")) + f = StringIO() + benzene.writepdb(f) + f.seek(0) + for i in range(2): + line = f.readline() + + refline = "ATOM 1 C LIG 0 1.194 -0.689 0.000 1.00 0.00 C \n" + assert line == refline + + +def test_mol2_format(xyz_folder, tmp_path): + """ + Test that a Mol2 file can be written and read for a molecule with double and aromatic bonds + """ + mol = Molecule(os.path.join(xyz_folder, "reactant2.xyz")) + mol.guess_bonds() + + f = StringIO() + mol.writemol2(f) + f.seek(0) + + mol2 = Molecule() + mol2.readmol2(f) + + assert len(mol) == len(mol2) + assert len(mol.bonds) == len(mol2.bonds) + + def test_read_multiple_molecules_from_coskf(coskf_folder): """ Test for COSKF reading