Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix bug in GRO reader #775

Merged
merged 10 commits into from
Jan 15, 2024
30 changes: 17 additions & 13 deletions gmso/formats/gro.py
Original file line number Diff line number Diff line change
Expand Up @@ -42,8 +42,8 @@ def read_gro(filename):
Gro files do not specify connections between atoms, the returned topology
will not have connections between sites either.

Currently this implementation does not support a gro file with more than 1
frame.
Currently this implementation does not support parsing velocities from a gro file or gro file
with more than 1 frame.

All residues and resid information from the gro file are currently lost
when converting to `topology`.
Expand All @@ -57,32 +57,36 @@ def read_gro(filename):
coords = u.nm * np.zeros(shape=(n_atoms, 3))
for row, _ in enumerate(coords):
line = gro_file.readline()
content = line.split()
if not line:
msg = (
"Incorrect number of lines in .gro file. Based on the "
"number in the second line of the file, {} rows of"
"atoms were expected, but at least one fewer was found."
)
raise ValueError(msg.format(n_atoms))
res_id = int(line[:5].strip())
res_name = line[5:10].strip()
atom_name = line[10:15].strip()
atom_id = line[15:20].strip()
Dismissed Show dismissed Hide dismissed

res = content[0]
atom_name = content[1]
atom_id = content[2]
positions = line[20:].split()
daico007 marked this conversation as resolved.
Show resolved Hide resolved
coords[row] = u.nm * np.array(
[
float(content[3]),
float(content[4]),
float(content[5]),
float(positions[0]),
float(positions[1]),
float(positions[2]),
]
)
site = Atom(name=atom_name, position=coords[row])

r = re.compile("([0-9]+)([a-zA-Z]+)")
m = r.match(res)
site.molecule = (m.group(2), int(m.group(1)))
site.residue = (m.group(2), int(m.group(1)))
site.molecule = (res_name, res_id)
site.residue = (res_name, res_id)
top.add_site(site, update_types=False)

if len(positions) == 6:
warnings.warn(
"Velocity information presents but will not be parsed."
)
top.update_topology()

# Box information
Expand Down
2 changes: 1 addition & 1 deletion gmso/tests/test_gro.py
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ def test_read_gro(self):
)

def test_wrong_n_atoms(self):
with pytest.raises(IndexError):
with pytest.raises(ValueError):
Topology.load(get_fn("too_few_atoms.gro"))
with pytest.raises(ValueError):
Topology.load(get_fn("too_many_atoms.gro"))
Expand Down
Loading