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.

Sign up for GitHub

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

Jump to bottom

Fixes AtomGroup.center when using compounds + unwrapping #2992

Merged
merged 1 commit into from
Oct 18, 2020
Merged

Fixes AtomGroup.center when using compounds + unwrapping #2992

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 3 additions & 1 deletion package/CHANGELOG
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,11 +15,13 @@ The rules for this file:
------------------------------------------------------------------------------
??/??/?? tylerjereddy, richardjgowers, IAlibay, hmacdope, orbeckst, cbouy,
lilyminium, daveminh, jbarnoud, yuxuanzhuang, VOD555, ianmkenney,
calcraven，xiki-tempula, mieczyslaw
calcraven，xiki-tempula, mieczyslaw, manuel.nuno.melo

* 2.0.0

Fixes
* AtomGroup.center now works correctly for compounds + unwrapping
(Issue #2984)
* Avoid using deprecated array indexing in topology attributes
(Issue #2990, PR #2991)
* ParmedParser no longer guesses elements if they are not recognised, instead
Expand Down
10 changes: 8 additions & 2 deletions package/MDAnalysis/core/groups.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -830,12 +830,18 @@ def center(self, weights, pbc=False, compound='group', unwrap=False):

# Sort positions and weights by compound index and promote to dtype if
# required:
sort_indices = np.argsort(compound_indices)

# are we already sorted? argsorting and fancy-indexing can be expensive
if np.any(np.diff(compound_indices) < 0):
sort_indices = np.argsort(compound_indices)
else:
sort_indices = slice(None)
compound_indices = compound_indices[sort_indices]

# Unwrap Atoms
if unwrap:
coords = atoms.unwrap(compound=comp, reference=None, inplace=False)
coords = atoms.unwrap(compound=comp, reference=None,
inplace=False)[sort_indices]
else:
coords = atoms.positions[sort_indices]
if weights is None:
Expand Down
67 changes: 45 additions & 22 deletions testsuite/MDAnalysisTests/core/test_atomgroup.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -991,6 +991,12 @@ def ag(self):
group.wrap(inplace=True)
return group

@pytest.fixture()
def unordered_ag(self, ag):
ndx = np.arange(len(ag))
np.random.shuffle(ndx)
return ag[ndx]

@pytest.fixture()
def ref_noUnwrap_residues(self):
return {
Expand All @@ -1013,12 +1019,12 @@ def ref_Unwrap_residues(self):
'COG': np.array([[21.356, 41.685, 40.501],
[44.577, 43.312, 79.039],
[ 2.204, 27.722, 54.023]], dtype=np.float32),
'COM': np.array([[20.815, 42.013, 39.802],
[44.918, 43.282, 79.325],
[2.045, 28.243, 54.127]], dtype=np.float32),
'MOI': np.array([[16747.486, -1330.489, 2938.243],
[-1330.489, 19315.253, 3306.212],
[ 2938.243, 3306.212, 8990.481]]),
'COM': np.array([[21.286, 41.664, 40.465],
[44.528, 43.426, 78.671],
[ 2.111, 27.871, 53.767]], dtype=np.float32),
'MOI': np.array([[16687.941, -1330.617, 2925.883],
[-1330.617, 19256.178, 3354.832],
[ 2925.883, 3354.832, 8989.946]]),
'Asph': 0.2969491080,
}

Expand Down Expand Up @@ -1058,6 +1064,12 @@ def test_UnWrapFlag_residues(self, ag, ref_Unwrap_residues):
assert_almost_equal(ag.moment_of_inertia(unwrap=True, compound='residues'), ref_Unwrap_residues['MOI'], self.prec)
assert_almost_equal(ag.asphericity(unwrap=True, compound='residues'), ref_Unwrap_residues['Asph'], self.prec)

def test_UnWrapFlag_residues_unordered(self, unordered_ag, ref_Unwrap_residues):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This could be parametrized so each AG method is tested separately.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You suggest also parameterizing the method? Like this?

@pytest.mark.parametrize('method_name', ('center_of_geometry', 'center_of_mass',
                                         'moment_of_inertia', 'asphericity'))
def test_UnWrapFlag_residues_unordered(self, unordered_ag, ref_Unwrap_residues, method_name):
    method = getattr(unordered_ag, method_name)
    # requires that the keys of ref_Unwrap_residues are changed to match the method names
    assert_almost_equal(method(unwrap=True, compound='residues'),
                        ref_Unwrap_residues[method_name], self.prec)

I'm all for it, but several more tests are begging for the same cleanup. It seems a targeted refactoring effort would make more sense. If you think that at least new tests should be parameterized as above, I'm cool with having it thus.

assert_almost_equal(unordered_ag.center_of_geometry(unwrap=True, compound='residues'), ref_Unwrap_residues['COG'], self.prec)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These lines could be split.

Copy link
Member Author

@mnmelo mnmelo Oct 17, 2020 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I really just kept them long to conform to how all other such tests are. I can split (shall I also split the analogous ones in neighboring tests?)

assert_almost_equal(unordered_ag.center_of_mass(unwrap=True, compound='residues'), ref_Unwrap_residues['COM'], self.prec)
assert_almost_equal(unordered_ag.moment_of_inertia(unwrap=True, compound='residues'), ref_Unwrap_residues['MOI'], self.prec)
assert_almost_equal(unordered_ag.asphericity(unwrap=True, compound='residues'), ref_Unwrap_residues['Asph'], self.prec)

def test_default(self, ref_noUnwrap):
u = UnWrapUniverse(is_triclinic=False)
group = u.atoms[31:39] # molecules 11
Expand Down Expand Up @@ -1281,22 +1293,27 @@ def test_center_of_mass_compounds(self, ag, name, compound):

@pytest.mark.parametrize('name, compound', (('resids', 'residues'),
('segids', 'segments')))
def test_center_of_geometry_compounds_pbc(self, ag, name, compound):
@pytest.mark.parametrize('unwrap', (True, False))
def test_center_of_geometry_compounds_pbc(self, ag, name, compound,
unwrap):
ag.dimensions = [50, 50, 50, 90, 90, 90]
ref = [a.center_of_geometry() for a in ag.groupby(name).values()]
ref = [a.center_of_geometry(unwrap=unwrap)
for a in ag.groupby(name).values()]
ref = distances.apply_PBC(np.asarray(ref, dtype=np.float32),
ag.dimensions)
cog = ag.center_of_geometry(pbc=True, compound=compound)
ag.dimensions)
cog = ag.center_of_geometry(pbc=True, compound=compound, unwrap=unwrap)
assert_almost_equal(cog, ref, decimal=5)

@pytest.mark.parametrize('name, compound', (('resids', 'residues'),
('segids', 'segments')))
def test_center_of_mass_compounds_pbc(self, ag, name, compound):
@pytest.mark.parametrize('unwrap', (True, False))
def test_center_of_mass_compounds_pbc(self, ag, name, compound, unwrap):
ag.dimensions = [50, 50, 50, 90, 90, 90]
ref = [a.center_of_mass() for a in ag.groupby(name).values()]
ref = [a.center_of_mass(unwrap=unwrap)
for a in ag.groupby(name).values()]
ref = distances.apply_PBC(np.asarray(ref, dtype=np.float32),
ag.dimensions)
com = ag.center_of_mass(pbc=True, compound=compound)
ag.dimensions)
com = ag.center_of_mass(pbc=True, compound=compound, unwrap=unwrap)
assert_almost_equal(com, ref, decimal=5)

@pytest.mark.parametrize('name, compound', (('molnums', 'molecules'),
Expand All @@ -1317,24 +1334,30 @@ def test_center_of_mass_compounds_special(self, ag_molfrg,

@pytest.mark.parametrize('name, compound', (('molnums', 'molecules'),
('fragindices', 'fragments')))
@pytest.mark.parametrize('unwrap', (True, False))
def test_center_of_geometry_compounds_special_pbc(self, ag_molfrg,
name, compound):
name, compound, unwrap):
ag_molfrg.dimensions = [50, 50, 50, 90, 90, 90]
ref = [a.center_of_geometry() for a in ag_molfrg.groupby(name).values()]
ref = [a.center_of_geometry(unwrap=unwrap)
for a in ag_molfrg.groupby(name).values()]
ref = distances.apply_PBC(np.asarray(ref, dtype=np.float32),
ag_molfrg.dimensions)
cog = ag_molfrg.center_of_geometry(pbc=True, compound=compound)
ag_molfrg.dimensions)
cog = ag_molfrg.center_of_geometry(pbc=True, compound=compound,
unwrap=unwrap)
assert_almost_equal(cog, ref, decimal=5)

@pytest.mark.parametrize('name, compound', (('molnums', 'molecules'),
('fragindices', 'fragments')))
@pytest.mark.parametrize('unwrap', (True, False))
def test_center_of_mass_compounds_special_pbc(self, ag_molfrg,
name, compound):
name, compound, unwrap):
ag_molfrg.dimensions = [50, 50, 50, 90, 90, 90]
ref = [a.center_of_mass() for a in ag_molfrg.groupby(name).values()]
ref = [a.center_of_mass(unwrap=unwrap)
for a in ag_molfrg.groupby(name).values()]
ref = distances.apply_PBC(np.asarray(ref, dtype=np.float32),
ag_molfrg.dimensions)
com = ag_molfrg.center_of_mass(pbc=True, compound=compound)
ag_molfrg.dimensions)
com = ag_molfrg.center_of_mass(pbc=True, compound=compound,
unwrap=unwrap)
assert_almost_equal(com, ref, decimal=5)

def test_center_wrong_compound(self, ag):
Expand Down