Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
17 changes: 14 additions & 3 deletions dpdata/md/rdf.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@
import numpy as np


def rdf(sys, sel_type=[None, None], max_r=5, nbins=100):
def rdf(sys, sel_type=None, max_r=5, nbins=100):
"""Compute the rdf of a system.

Parameters
Expand Down Expand Up @@ -40,7 +40,12 @@ def rdf(sys, sel_type=[None, None], max_r=5, nbins=100):
)


def compute_rdf(box, posis, atype, sel_type=[None, None], max_r=5, nbins=100):
def compute_rdf(box, posis, atype, sel_type=None, max_r=5, nbins=100):
"""Compute an RDF without mutating the caller's selection list."""
if sel_type is None:
sel_type = [None, None]
else:
sel_type = list(sel_type)
nframes = box.shape[0]
xx = None
all_rdf = []
Expand All @@ -58,7 +63,13 @@ def compute_rdf(box, posis, atype, sel_type=[None, None], max_r=5, nbins=100):
return xx, all_rdf, all_cod


def _compute_rdf_1frame(box, posis, atype, sel_type=[None, None], max_r=5, nbins=100):
def _compute_rdf_1frame(box, posis, atype, sel_type=None, max_r=5, nbins=100):
if sel_type is None:
sel_type = [None, None]
else:
# Normalise into a fresh list because the two ``None`` replacements
# below must not leak into a caller-owned list or a function default.
sel_type = list(sel_type)
all_types = list(set(list(np.sort(atype, kind="stable"))))
if sel_type[0] is None:
sel_type[0] = all_types
Expand Down
34 changes: 34 additions & 0 deletions tests/test_rdf.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@
from __future__ import annotations

import unittest

import numpy as np

from dpdata.md.rdf import compute_rdf

try:
import ase # noqa: F401
except ModuleNotFoundError:
skip_ase = True
else:
skip_ase = False


@unittest.skipIf(skip_ase, "RDF calculation requires ASE")
class TestRDFSelectionOwnership(unittest.TestCase):
def test_default_selection_is_recomputed_per_system(self):
box = np.eye(3).reshape(1, 3, 3) * 10
positions = np.array([[[0.0, 0.0, 0.0], [1.0, 0.0, 0.0]]])
rdf = compute_rdf

rdf(box, positions, np.array([0, 1]))
self.assertEqual(rdf.__defaults__[0], None)
_, values, _ = rdf(box, positions, np.array([2, 2]))
self.assertTrue(np.isfinite(values).all())

def test_caller_selection_is_not_mutated(self):
box = np.eye(3).reshape(1, 3, 3) * 10
positions = np.array([[[0.0, 0.0, 0.0], [1.0, 0.0, 0.0]]])
selection = [None, None]
compute_rdf(box, positions, np.array([0, 1]), selection)
self.assertEqual(selection, [None, None])
Loading