From e47ed4dadeb99904a990260c5fb5ef15e4c7d9e9 Mon Sep 17 00:00:00 2001 From: ayushsuhane Date: Thu, 14 Jun 2018 18:25:39 -0700 Subject: [PATCH 01/17] Clean initial Capped distance feature --- package/MDAnalysis/lib/distances.py | 112 ++++++++++++++++++++++++++++ 1 file changed, 112 insertions(+) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index 11ba3c1e112..b16d1a66b6c 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -389,6 +389,118 @@ def self_distance_array(reference, box=None, result=None, backend="serial"): return distances +def capped_distance(reference, configuration, max_cutoff, min_cutoff = 0, box = None, method = None): + """Returns the indices of reference and configuration atoms within a fixed distance + + If a box is supplied, then a minimum image convention is used to evaluate + the distances. + + An optional keyword for the method is provided to allow usage of other + implemented methods. Currently pkdtree and bruteforce can be used. Depending + on the number of particles and distribution, the times may vary significantly. + + Parameters + ----------- + reference : numpy.array of type numpy.float32 + The shape of the reference array should be ``(N, 3)`` + or a (3, ) coordinate. + configuration : numpy.array of type numpy.float32 + All the attributes are similar to reference array. Essentially + second group of atoms. + max_cutoff : float32 + The maximum cutoff distance. Only the pair of indices which + have smaller distance than this cutoff distance will be returned. + min_cutoff : float32 + The lower bound of cutoff. The function will not return the + pairs of indices whose distances are less than this parameter. + box : (optional) array or None, default ``None`` + The dimensions, if provided, must be provided in the same + The unitcell dimesions for this system format as returned + by :attr:`MDAnalysis.coordinates.base.Timestep.dimensions`: + ``[lx,ly, lz, alpha, beta, gamma]``. Minimum image convention + is applied if the box is provided. + method : 'bruteforce' or 'pkdtree' or 'None' + Based on benchmarks, checks the suitable method + and search the indices of the selections in relatively + less amount of time. 'bruteforce' can be costly with + large number of coordinates. Another option is + periodic KDTree and can be advantageous with increase in particles. + + Returns + ------- + indices : list + Combination of coordinate indices, one from each + from reference and configuration state such that + distance between them is smaller than cutoff distance + provided as ``cutoff`` in arguments. + + Note + ----- + Currently only supports brute force. Similar methods + can be defined and can be directly linked + to this function. + """ + + if box is not None: + boxtype = _box_check(box) + # Convert [A,B,C,alpha,beta,gamma] to [[A],[B],[C]] + if (boxtype == 'tri_box'): + box = triclinic_vectors(box) + if (boxtype == 'tri_vecs_bad'): + box = triclinic_vectors(triclinic_box(box[0], box[1], box[2])) + if method is None: + method = _determine_method(reference, configuration, max_cutoff, box = box) + + if method == 'bruteforce': + pairs = _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff = min_cutoff, box = box) + #if method = 'pkdt': + # pairs = _pkdt_capped(reference, configuration, cutoff, box = None) + return pairs + + +def _determine_method(reference, configuration, cutoff, box): + """ + Switch between different methods based on the the optimized time. + All the rules to select the method based on the input can be + incorporated here. + + Returns + ------- + String of the method. supports ``bruteforce`` or ``pkdt`` as of now + + Note + ---- + Just a basic rule based on the number of particles + """ + #if (len(reference) > 1000) or (len(configuration) > 1000): + # return 'pkdt' + return 'bruteforce' + + +def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff = 0, box = None): + """ + Using naive distanc calulations, returns a list + containing the indices with one from each + reference and configuration list, such that the distance between + them is less than the specified cutoff distance + """ + pairs = [] + if reference.shape == (3,): + reference = reference[None, :] + if configuration.shape == (3,): + configuration = configuration[None, :] + + _check_array(reference, 'reference') + _check_array(configuration, 'configuration') + + for i, coords in enumerate(reference): + dist = distance_array(coords[None, :], configuration, box=box)[0] + idx = np.where((dist <= max_cutoff) & (dist > min_cutoff))[0] + for j in idx: + pairs.append((i,j)) + + return pairs + def transform_RtoS(inputcoords, box, backend="serial"): """Transform an array of coordinates from real space to S space (aka lambda space) From 23e9868801d9e572039d32c7db7f9659029a8649 Mon Sep 17 00:00:00 2001 From: ayush Date: Sun, 24 Jun 2018 15:36:27 -0700 Subject: [PATCH 02/17] Modification as per the review --- package/MDAnalysis/lib/distances.py | 142 ++++++++++-------- .../MDAnalysisTests/lib/test_distances.py | 52 +++++++ 2 files changed, 134 insertions(+), 60 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index b16d1a66b6c..bc37264192d 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -389,58 +389,70 @@ def self_distance_array(reference, box=None, result=None, backend="serial"): return distances -def capped_distance(reference, configuration, max_cutoff, min_cutoff = 0, box = None, method = None): - """Returns the indices of reference and configuration atoms within a fixed distance +def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=None, method=None): + """Calculates the pairs and distance within a specified distance - If a box is supplied, then a minimum image convention is used to evaluate - the distances. + If a *box* is supplied, then a minimum image convention is used + to evaluate the distances. + + An automatic guessing of optimized method to calculate the distances is + included in the function. An optional keyword for the method is also + provided. Users can override the method with this functionality. + Currently pkdtree and bruteforce are implemented. - An optional keyword for the method is provided to allow usage of other - implemented methods. Currently pkdtree and bruteforce can be used. Depending - on the number of particles and distribution, the times may vary significantly. Parameters ----------- - reference : numpy.array of type numpy.float32 - The shape of the reference array should be ``(N, 3)`` - or a (3, ) coordinate. - configuration : numpy.array of type numpy.float32 - All the attributes are similar to reference array. Essentially - second group of atoms. - max_cutoff : float32 - The maximum cutoff distance. Only the pair of indices which - have smaller distance than this cutoff distance will be returned. - min_cutoff : float32 - The lower bound of cutoff. The function will not return the - pairs of indices whose distances are less than this parameter. - box : (optional) array or None, default ``None`` + reference : array + The shape of the reference array should be ``(N, 3)`` + or a (3, ) coordinate. + configuration : array + Similar to reference array and second group of coordinates + max_cutoff : float + Maximum cutoff distance between the reference and configuration + min_cutoff : (optional) float + Minimum cutoff distance between reference and configuration [None] + box : (optional) array or None The dimensions, if provided, must be provided in the same The unitcell dimesions for this system format as returned by :attr:`MDAnalysis.coordinates.base.Timestep.dimensions`: ``[lx,ly, lz, alpha, beta, gamma]``. Minimum image convention - is applied if the box is provided. - method : 'bruteforce' or 'pkdtree' or 'None' - Based on benchmarks, checks the suitable method - and search the indices of the selections in relatively - less amount of time. 'bruteforce' can be costly with - large number of coordinates. Another option is - periodic KDTree and can be advantageous with increase in particles. - + is applied if the box is provided [None] + method : (optional) 'bruteforce' or 'pkdtree' or 'None' + Keyword to override the automatic guessing of method built-in + in the function [None] + Returns ------- - indices : list - Combination of coordinate indices, one from each - from reference and configuration state such that - distance between them is smaller than cutoff distance - provided as ``cutoff`` in arguments. - + pairs : array + Pair of indices, one from each reference and configuration such that + distance between them is within the ``max_cutoff`` and ``min_cutoff`` + pairs[i,j] contains the indices i from reference coordinates, and + j from configuration + distances : array + Distances corresponding to each pair of indices. + d[k] corresponding to the pairs[i,j] gives the distance between + i-th and j-th coordinate in reference and configuration respectively + + .. code-block:: python + + pairs, distances = capped_distances(reference, coordinates, max_cutoff) + for indx, [a,b] in enumerate(pairs): + coord1 = reference[a] + coord2 = configuration[b] + distance = distances[indx] + Note ----- - Currently only supports brute force. Similar methods - can be defined and can be directly linked - to this function. + Currently only supports brute force. + + .. SeeAlso:: :func:'MDAnalysis.lib.distances.distance_array' + + Similar methods can be defined and can be directly linked to this function. + + """ - + if box is not None: boxtype = _box_check(box) # Convert [A,B,C,alpha,beta,gamma] to [[A],[B],[C]] @@ -449,45 +461,52 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff = 0, box = if (boxtype == 'tri_vecs_bad'): box = triclinic_vectors(triclinic_box(box[0], box[1], box[2])) if method is None: - method = _determine_method(reference, configuration, max_cutoff, box = box) + method = _determine_method(reference, configuration, + max_cutoff, box=box) if method == 'bruteforce': - pairs = _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff = min_cutoff, box = box) - #if method = 'pkdt': - # pairs = _pkdt_capped(reference, configuration, cutoff, box = None) - return pairs - + pairs, dist = _bruteforce_capped(reference, configuration, + max_cutoff, min_cutoff=min_cutoff, + box=box) + # if method = 'pkdt': + # pairs = _pkdt_capped(reference, configuration, cutoff, box = None) + return np.array(pairs), np.array(dist) + def _determine_method(reference, configuration, cutoff, box): """ Switch between different methods based on the the optimized time. - All the rules to select the method based on the input can be + All the rules to select the method based on the input can be incorporated here. Returns ------- - String of the method. supports ``bruteforce`` or ``pkdt`` as of now + String of the method. Supports ``bruteforce`` or ``pkdt`` as of now - Note + Note ---- Just a basic rule based on the number of particles """ - #if (len(reference) > 1000) or (len(configuration) > 1000): - # return 'pkdt' + # if (len(reference) > 1000) or (len(configuration) > 1000): + # return 'pkdt' return 'bruteforce' -def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff = 0, box = None): - """ - Using naive distanc calulations, returns a list - containing the indices with one from each +def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=None, box=None): + """ + Using naive distance calulations, returns a list + containing the indices with one from each reference and configuration list, such that the distance between them is less than the specified cutoff distance """ - pairs = [] - if reference.shape == (3,): + pairs, distance = [], [] + + reference = np.asarray(reference, dtype=np.float32) + configuration = np.asarray(configuration, dtype=np.float32) + + if reference.shape == (3, ): reference = reference[None, :] - if configuration.shape == (3,): + if configuration.shape == (3, ): configuration = configuration[None, :] _check_array(reference, 'reference') @@ -495,11 +514,14 @@ def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff = 0, bo for i, coords in enumerate(reference): dist = distance_array(coords[None, :], configuration, box=box)[0] - idx = np.where((dist <= max_cutoff) & (dist > min_cutoff))[0] + if min_cutoff is not None: + idx = np.where((dist <= max_cutoff) & (dist > min_cutoff))[0] + else: + idx = np.where((dist <= max_cutoff))[0] for j in idx: - pairs.append((i,j)) - - return pairs + pairs.append((i, j)) + distance.append(dist[j]) + return pairs, distance def transform_RtoS(inputcoords, box, backend="serial"): diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index 71e74231239..5061bec2475 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -27,6 +27,8 @@ import MDAnalysis as mda +from MDAnalysis.lib.mdamath import triclinic_vectors + def test_transform_StoR_pass(): box = np.array([10, 7, 3, 45, 60, 90], dtype=np.float32) s = np.array([[0.5, -0.1, 0.5]], dtype=np.float32) @@ -37,9 +39,59 @@ def test_transform_StoR_pass(): assert_equal(original_r, test_r) + def test_transform_StoR_fail(): box = np.array([10, 7, 3, 45, 60, 90], dtype=np.float32) s = np.array([[0.5, -0.1, 0.5]]) with pytest.raises(TypeError, match='S must be of type float32'): r = mda.lib.distances.transform_StoR(s, box) + + +def test_capped_distance_noresults(): + point1 = np.array([0.1, 0.1, 0.1], dtype=np.float32) + point2 = np.array([0.95, 0.1, 0.1], dtype=np.float32) + + pairs, distances = mda.lib.distances.capped_distance(point1, + point2, + max_cutoff=0.2) + + assert_equal(len(pairs), 0) + + +boxes_1 = (np.array([1, 2, 3, 90, 90, 90], dtype=np.float32), # ortho + np.array([1, 2, 3, 30, 45, 60], dtype=np.float32), # tri_box + triclinic_vectors( # tri_vecs + np.array([1, 2, 3, 90, 90, 45], dtype=np.float32)), + np.array([[0.5, 0.9, 1.9], # tri_vecs_bad + [2.0, 0.4, 0.1], + [0.0, 0.6, 0.5]], dtype=np.float32)) + +query_1 = (np.array([0.1, 0.1, 0.1], dtype=np.float32), + np.array([[0.1, 0.1, 0.1], + [0.2, 0.1, 0.1]], dtype=np.float32)) + +points = (np.random.uniform(low=0, high=1.0, + size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) + +np.random.seed(90003) + + +@pytest.mark.parametrize('box, query', zip(boxes_1, query_1)) +def test_capped_distance_checkbrute(box, query): + max_cutoff = 0.3 + # capped distance should be able to handle array of vectors + # as well as single vectors. + pairs, dist = mda.lib.distances.capped_distance(query, + points, + max_cutoff, + box=box) + if(query.shape[0] == 3): + distances = mda.lib.distances.distance_array(query[None, :], + points, box=box) + else: + distances = mda.lib.distances.distance_array(query, + points, box=box) + indices = np.where(distances < max_cutoff) + + assert_equal(pairs[:, 1], indices[1]) From dd12343150a16a53dcf7e05a6894a1c1965dcb7c Mon Sep 17 00:00:00 2001 From: ayush Date: Sun, 24 Jun 2018 16:29:39 -0700 Subject: [PATCH 03/17] modified zip with itertools.product --- testsuite/MDAnalysisTests/lib/test_distances.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index 5061bec2475..cc0f01c2017 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -25,6 +25,8 @@ import numpy as np from numpy.testing import assert_equal +import itertools + import MDAnalysis as mda from MDAnalysis.lib.mdamath import triclinic_vectors @@ -77,7 +79,7 @@ def test_capped_distance_noresults(): np.random.seed(90003) -@pytest.mark.parametrize('box, query', zip(boxes_1, query_1)) +@pytest.mark.parametrize('box, query', itertools.product(boxes_1, query_1)) def test_capped_distance_checkbrute(box, query): max_cutoff = 0.3 # capped distance should be able to handle array of vectors From 9833a6909d674bcf7dc35415cc2133e25d6f9c07 Mon Sep 17 00:00:00 2001 From: ayush Date: Sun, 24 Jun 2018 17:33:28 -0700 Subject: [PATCH 04/17] random seed --- testsuite/MDAnalysisTests/lib/test_distances.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index cc0f01c2017..0eebcdcbb46 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -73,10 +73,11 @@ def test_capped_distance_noresults(): np.array([[0.1, 0.1, 0.1], [0.2, 0.1, 0.1]], dtype=np.float32)) +np.random.seed(90003) points = (np.random.uniform(low=0, high=1.0, size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) -np.random.seed(90003) + @pytest.mark.parametrize('box, query', itertools.product(boxes_1, query_1)) From e1dce21d0cfb28446a2bda05a5fc23ced34c4e0e Mon Sep 17 00:00:00 2001 From: ayush Date: Sun, 24 Jun 2018 19:56:46 -0700 Subject: [PATCH 05/17] Capped PKDTree --- package/MDAnalysis/lib/distances.py | 65 +++++++++++++++++++++++++++-- 1 file changed, 62 insertions(+), 3 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index bc37264192d..420af741ee8 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -72,6 +72,7 @@ from .mdamath import triclinic_vectors, triclinic_box + # hack to select backend with backend= kwarg. Note that # the cython parallel code (prange) in parallel.distances is # independent from the OpenMP code @@ -389,6 +390,7 @@ def self_distance_array(reference, box=None, result=None, backend="serial"): return distances + def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=None, method=None): """Calculates the pairs and distance within a specified distance @@ -444,9 +446,10 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N Note ----- - Currently only supports brute force. + Currently only supports brute force and Periodic KDtree .. SeeAlso:: :func:'MDAnalysis.lib.distances.distance_array' + .. SeeAlso:: :func:'MDAnalysis.lib.pkdtree.PeriodicKDTree' Similar methods can be defined and can be directly linked to this function. @@ -468,8 +471,10 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N pairs, dist = _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=min_cutoff, box=box) - # if method = 'pkdt': - # pairs = _pkdt_capped(reference, configuration, cutoff, box = None) + if method == 'pkdt': + pairs, dist = _pkdt_capped(reference, configuration, + max_cutoff, min_cutoff=min_cutoff, + box=box) return np.array(pairs), np.array(dist) @@ -523,6 +528,60 @@ def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=None, bo distance.append(dist[j]) return pairs, distance +def _pkdtree_capped(reference, configuration, max_cutoff, min_cutoff=None, box=None): + """ Capped Distance evaluations using KDtree. + + Uses minimum image convention if *box* is specified + + Returns: + -------- + pairs : list + List of atom indices which are within the specified cutoff distance. + pairs `(i, j)` corresponds to i-th particle in reference and + j-th particle in configuration + distance : list + Distance between two atoms corresponding to the (i, j) indices + in pairs. + + """ + from .pkdtree import PeriodicKDTree + from Bio.KDTree import KDTree + + pairs, distances = [], [] + + reference = np.asarray(reference, dtype=np.float32) + configuration = np.asarray(configuration, dtype=np.float32) + + if reference.shape == (3, ): + reference = reference[None, :] + if configuration.shape == (3, ): + configuration = configuration[None, :] + + _check_array(reference, 'reference') + _check_array(configuration, 'configuration') + + # Build The KDTree + if box is not None: + kdtree = PeriodicKDTree(box, bucket_size=10) + else : + kdtree = KDTree(dim=3, bucket_size=10) + + kdtree.set_coords(configuration) + # Search for every query point + for idx,centers in enumerate(reference): + kdtree.search(centers, max_cutoff) + indices = kdtree.get_indices() + dist = distance_array(centers.reshape((1,3)), configuration[indices],box=box)[0] + if min_cutoff is not None: + mask = np.where(dist > min_cutoff)[0] + dist = dist[mask] + indices = [indices[mask[i]] for i in range(len(mask))] + if len(indices) != 0: + for num, j in enumerate(indices): + pairs.append((idx,j)) + distances.append(dist[num]) + return pairs, distances + def transform_RtoS(inputcoords, box, backend="serial"): """Transform an array of coordinates from real space to S space (aka lambda space) From 8cc25fe2df0e12d9910440f7cb7664ad52f6c6b7 Mon Sep 17 00:00:00 2001 From: ayush Date: Sun, 24 Jun 2018 20:33:18 -0700 Subject: [PATCH 06/17] PEP8 and added test --- package/MDAnalysis/lib/distances.py | 24 ++++++++++--------- .../MDAnalysisTests/lib/test_distances.py | 12 ++++++---- 2 files changed, 20 insertions(+), 16 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index 420af741ee8..8c772b55d6a 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -471,10 +471,10 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N pairs, dist = _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=min_cutoff, box=box) - if method == 'pkdt': - pairs, dist = _pkdt_capped(reference, configuration, - max_cutoff, min_cutoff=min_cutoff, - box=box) + if method == 'pkdtree': + pairs, dist = _pkdtree_capped(reference, configuration, + max_cutoff, min_cutoff=min_cutoff, + box=box) return np.array(pairs), np.array(dist) @@ -493,12 +493,12 @@ def _determine_method(reference, configuration, cutoff, box): Just a basic rule based on the number of particles """ # if (len(reference) > 1000) or (len(configuration) > 1000): - # return 'pkdt' + # return 'pkdtree' return 'bruteforce' def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=None, box=None): - """ + """ Using naive distance calulations, returns a list containing the indices with one from each reference and configuration list, such that the distance between @@ -528,6 +528,7 @@ def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=None, bo distance.append(dist[j]) return pairs, distance + def _pkdtree_capped(reference, configuration, max_cutoff, min_cutoff=None, box=None): """ Capped Distance evaluations using KDtree. @@ -537,7 +538,7 @@ def _pkdtree_capped(reference, configuration, max_cutoff, min_cutoff=None, box=N -------- pairs : list List of atom indices which are within the specified cutoff distance. - pairs `(i, j)` corresponds to i-th particle in reference and + pairs `(i, j)` corresponds to i-th particle in reference and j-th particle in configuration distance : list Distance between two atoms corresponding to the (i, j) indices @@ -563,22 +564,23 @@ def _pkdtree_capped(reference, configuration, max_cutoff, min_cutoff=None, box=N # Build The KDTree if box is not None: kdtree = PeriodicKDTree(box, bucket_size=10) - else : + else: kdtree = KDTree(dim=3, bucket_size=10) kdtree.set_coords(configuration) # Search for every query point - for idx,centers in enumerate(reference): + for idx, centers in enumerate(reference): kdtree.search(centers, max_cutoff) indices = kdtree.get_indices() - dist = distance_array(centers.reshape((1,3)), configuration[indices],box=box)[0] + dist = distance_array(centers.reshape((1, 3)), + configuration[indices], box=box)[0] if min_cutoff is not None: mask = np.where(dist > min_cutoff)[0] dist = dist[mask] indices = [indices[mask[i]] for i in range(len(mask))] if len(indices) != 0: for num, j in enumerate(indices): - pairs.append((idx,j)) + pairs.append((idx, j)) distances.append(dist[num]) return pairs, distances diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index 0eebcdcbb46..1813ec9c275 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -73,22 +73,24 @@ def test_capped_distance_noresults(): np.array([[0.1, 0.1, 0.1], [0.2, 0.1, 0.1]], dtype=np.float32)) +method_1 = ('bruteforce', 'pkdtree') + np.random.seed(90003) points = (np.random.uniform(low=0, high=1.0, size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) - - -@pytest.mark.parametrize('box, query', itertools.product(boxes_1, query_1)) -def test_capped_distance_checkbrute(box, query): +@pytest.mark.parametrize('box, query , method', + itertools.product(boxes_1, query_1, method_1)) +def test_capped_distance_checkbrute(box, query, method): max_cutoff = 0.3 # capped distance should be able to handle array of vectors # as well as single vectors. pairs, dist = mda.lib.distances.capped_distance(query, points, max_cutoff, - box=box) + box=box, + method=method) if(query.shape[0] == 3): distances = mda.lib.distances.distance_array(query[None, :], points, box=box) From f02868cf7234761e3e31c757f48569b508e72289 Mon Sep 17 00:00:00 2001 From: ayush Date: Tue, 26 Jun 2018 14:07:42 -0700 Subject: [PATCH 07/17] Added rules for method selection and modifications --- package/MDAnalysis/lib/distances.py | 70 +++++++++++++++++------------ 1 file changed, 42 insertions(+), 28 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index 8c772b55d6a..341a1307d6b 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -1,4 +1,4 @@ -# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; -*- + # -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; -*- # vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 # # MDAnalysis --- https://www.mdanalysis.org @@ -406,10 +406,11 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N Parameters ----------- reference : array - The shape of the reference array should be ``(N, 3)`` - or a (3, ) coordinate. + reference coordinates array with shape ``reference.shape = (3,)`` + or ``reference.shape = (len(reference), 3)`` configuration : array - Similar to reference array and second group of coordinates + Configuration coordinate array with shape ``reference.shape = (3,)`` + or ``reference.shape = (len(reference), 3)`` max_cutoff : float Maximum cutoff distance between the reference and configuration min_cutoff : (optional) float @@ -451,9 +452,6 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N .. SeeAlso:: :func:'MDAnalysis.lib.distances.distance_array' .. SeeAlso:: :func:'MDAnalysis.lib.pkdtree.PeriodicKDTree' - Similar methods can be defined and can be directly linked to this function. - - """ if box is not None: @@ -463,22 +461,16 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N box = triclinic_vectors(box) if (boxtype == 'tri_vecs_bad'): box = triclinic_vectors(triclinic_box(box[0], box[1], box[2])) - if method is None: - method = _determine_method(reference, configuration, - max_cutoff, box=box) - - if method == 'bruteforce': - pairs, dist = _bruteforce_capped(reference, configuration, - max_cutoff, min_cutoff=min_cutoff, - box=box) - if method == 'pkdtree': - pairs, dist = _pkdtree_capped(reference, configuration, - max_cutoff, min_cutoff=min_cutoff, - box=box) + method = _determine_method(reference, configuration, + max_cutoff, min_cutoff=min_cutoff, + box=box, method=method) + pairs, dist = method(reference, configuration, max_cutoff, + min_cutoff=min_cutoff, box=box) + return np.array(pairs), np.array(dist) -def _determine_method(reference, configuration, cutoff, box): +def _determine_method(reference, configuration, max_cutoff, min_cutoff=None, box=None, method=None): """ Switch between different methods based on the the optimized time. All the rules to select the method based on the input can be @@ -486,22 +478,44 @@ def _determine_method(reference, configuration, cutoff, box): Returns ------- - String of the method. Supports ``bruteforce`` or ``pkdt`` as of now + Function object based on the rules and specified method + + Currently implemented methods are + bruteforce : returns ``_bruteforce_capped`` + PKDtree : return ``_pkdtree_capped` - Note - ---- - Just a basic rule based on the number of particles """ - # if (len(reference) > 1000) or (len(configuration) > 1000): - # return 'pkdtree' - return 'bruteforce' + methods = {'bruteforce': _bruteforce_capped, + 'pkdtree':_pkdtree_capped} + + if method is not None: + return methods[method] + + if len(reference) > 5000 and len(configuration) > 5000: + if box is None and reference.shape[0] != 3 and configuration.shape[0] != 3: + min_dim = np.array([reference.min(axis=0), + configuration.min(axis=0)]) + max_dim = np.array([reference.max(axis=0), + configuration.max(axis=0)]) + size = max_dim.max(axis=0) - min_dim.min(axis=0) + elif box is not None: + if box.shape[0] == 6: + size = box[:3] + else: + size = box.max(axis=0) - box.min(axis=0) + + if np.any(size < 10.0*max_cutoff) and len(reference) > 100000 and len(configuration) > 100000: + return methods['bruteforce'] + else: + return methods['pkdtree'] + return methods['bruteforce'] def _bruteforce_capped(reference, configuration, max_cutoff, min_cutoff=None, box=None): """ Using naive distance calulations, returns a list containing the indices with one from each - reference and configuration list, such that the distance between + reference and configuration arrays, such that the distance between them is less than the specified cutoff distance """ pairs, distance = [], [] From 61d9c6b07096ac57f5be7ef682d8a5c78d5d3cf5 Mon Sep 17 00:00:00 2001 From: ayush Date: Sat, 30 Jun 2018 02:15:23 -0700 Subject: [PATCH 08/17] Changes the return type --- package/MDAnalysis/lib/distances.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index 341a1307d6b..51335f96dbb 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -466,8 +466,8 @@ def capped_distance(reference, configuration, max_cutoff, min_cutoff=None, box=N box=box, method=method) pairs, dist = method(reference, configuration, max_cutoff, min_cutoff=min_cutoff, box=box) - - return np.array(pairs), np.array(dist) + + return np.asarray(pairs), np.asarray(dist) def _determine_method(reference, configuration, max_cutoff, min_cutoff=None, box=None, method=None): From 57363bb6bb92365602b9f41d6b168bc376300303 Mon Sep 17 00:00:00 2001 From: ayush Date: Sat, 30 Jun 2018 12:01:34 -0700 Subject: [PATCH 09/17] fixed the testcase --- testsuite/MDAnalysisTests/lib/test_distances.py | 8 -------- 1 file changed, 8 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index 5f5883283bb..c9916b117db 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -43,14 +43,6 @@ def test_transform_StoR_pass(coord_dtype): assert_equal(original_r, test_r) -def test_transform_StoR_fail(): - box = np.array([10, 7, 3, 45, 60, 90], dtype=np.float32) - s = np.array([[0.5, -0.1, 0.5]]) - - with pytest.raises(TypeError, match='S must be of type float32'): - r = mda.lib.distances.transform_StoR(s, box) - - def test_capped_distance_noresults(): point1 = np.array([0.1, 0.1, 0.1], dtype=np.float32) point2 = np.array([0.95, 0.1, 0.1], dtype=np.float32) From ae6c2d991bac3c43838bf6be12ff186cb04d640a Mon Sep 17 00:00:00 2001 From: ayush Date: Tue, 3 Jul 2018 11:58:36 -0700 Subject: [PATCH 10/17] Stacked decorators in tests --- package/MDAnalysis/lib/distances.py | 11 ++++++----- testsuite/MDAnalysisTests/lib/test_distances.py | 8 ++++++-- 2 files changed, 12 insertions(+), 7 deletions(-) diff --git a/package/MDAnalysis/lib/distances.py b/package/MDAnalysis/lib/distances.py index 93ae5c45a52..790e540ed08 100644 --- a/package/MDAnalysis/lib/distances.py +++ b/package/MDAnalysis/lib/distances.py @@ -1,4 +1,4 @@ - # -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; -*- +# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; -*- # vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 # # MDAnalysis --- https://www.mdanalysis.org @@ -485,15 +485,14 @@ def _determine_method(reference, configuration, max_cutoff, min_cutoff=None, box Returns ------- Function object based on the rules and specified method - Currently implemented methods are bruteforce : returns ``_bruteforce_capped`` PKDtree : return ``_pkdtree_capped` """ methods = {'bruteforce': _bruteforce_capped, - 'pkdtree':_pkdtree_capped} - + 'pkdtree': _pkdtree_capped} + if method is not None: return methods[method] @@ -510,7 +509,9 @@ def _determine_method(reference, configuration, max_cutoff, min_cutoff=None, box else: size = box.max(axis=0) - box.min(axis=0) - if np.any(size < 10.0*max_cutoff) and len(reference) > 100000 and len(configuration) > 100000: + if (np.any(size < 10.0*max_cutoff) and + len(reference) > 100000 and + len(configuration) > 100000): return methods['bruteforce'] else: return methods['pkdtree'] diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index c9916b117db..0b29dd7f8eb 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -73,9 +73,13 @@ def test_capped_distance_noresults(): size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) -@pytest.mark.parametrize('box, query , method', - itertools.product(boxes_1, query_1, method_1)) +@pytest.mark.parametrize('box', boxes_1) +@pytest.mark.parametrize('query', query_1) +@pytest.mark.parametrize('method', method_1) def test_capped_distance_checkbrute(box, query, method): + np.random.seed(90003) + points = (np.random.uniform(low=0, high=1.0, + size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) max_cutoff = 0.3 # capped distance should be able to handle array of vectors # as well as single vectors. From 32b2f06e8afb8aa868f51562802a3a8aac9da121 Mon Sep 17 00:00:00 2001 From: ayush Date: Wed, 4 Jul 2018 08:50:48 -0700 Subject: [PATCH 11/17] Removed the seed --- testsuite/MDAnalysisTests/lib/test_distances.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index 0b29dd7f8eb..d2a87cbdfd7 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -68,11 +68,6 @@ def test_capped_distance_noresults(): method_1 = ('bruteforce', 'pkdtree') -np.random.seed(90003) -points = (np.random.uniform(low=0, high=1.0, - size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) - - @pytest.mark.parametrize('box', boxes_1) @pytest.mark.parametrize('query', query_1) @pytest.mark.parametrize('method', method_1) From 10685787f641fbab6ece1349b85df3551aad4948 Mon Sep 17 00:00:00 2001 From: ayush Date: Tue, 10 Jul 2018 18:15:05 -0700 Subject: [PATCH 12/17] Modified Changelog --- package/CHANGELOG | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/package/CHANGELOG b/package/CHANGELOG index d54408d699a..0f72cd2cb7d 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -14,12 +14,14 @@ The rules for this file: ------------------------------------------------------------------------------ ??/??/18 tylerjereddy, richardjgowers, palnabarun, orbeckst, kain88-de, zemanj, - VOD555, davidercruz, jbarnoud + VOD555, davidercruz, jbarnoud, ayushsuhane * 0.18.1 Enhancements + * Added capped_distance function for automatic selection of + distance evaluation method (PR #1941) * Added a rotation coordinate transformation (PR #1937) * Added a box centering trajectory transformation (PR #1946) * Added a on-the-fly trajectory transformations API and a coordinate translation From 221868656d2acab93829ec8b30167a1fc71b8bfb Mon Sep 17 00:00:00 2001 From: Richard Gowers Date: Wed, 11 Jul 2018 15:05:31 -0500 Subject: [PATCH 13/17] Update CHANGELOG --- package/CHANGELOG | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/package/CHANGELOG b/package/CHANGELOG index 0f72cd2cb7d..207eeb6dec5 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -20,8 +20,8 @@ The rules for this file: Enhancements - * Added capped_distance function for automatic selection of - distance evaluation method (PR #1941) + * Added lib.distances.capped_distance function to quickly calculate all distances + up to a given maximum distance (PR #1941) * Added a rotation coordinate transformation (PR #1937) * Added a box centering trajectory transformation (PR #1946) * Added a on-the-fly trajectory transformations API and a coordinate translation From cea2c379e5f82746d1c9390a2e32906d34ed1525 Mon Sep 17 00:00:00 2001 From: ayush Date: Wed, 11 Jul 2018 19:35:12 -0700 Subject: [PATCH 14/17] Added tests for minimum cutoff --- .../MDAnalysisTests/lib/test_distances.py | 20 +++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index d2a87cbdfd7..e9ff40f6777 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -68,10 +68,14 @@ def test_capped_distance_noresults(): method_1 = ('bruteforce', 'pkdtree') +min_cutoff_1 = (None, 0.1) + + @pytest.mark.parametrize('box', boxes_1) @pytest.mark.parametrize('query', query_1) @pytest.mark.parametrize('method', method_1) -def test_capped_distance_checkbrute(box, query, method): +@pytest.mark.parametrize('min_cutoff', min_cutoff_1) +def test_capped_distance_checkbrute(box, query, method, min_cutoff): np.random.seed(90003) points = (np.random.uniform(low=0, high=1.0, size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) @@ -81,15 +85,19 @@ def test_capped_distance_checkbrute(box, query, method): pairs, dist = mda.lib.distances.capped_distance(query, points, max_cutoff, + min_cutoff=min_cutoff, box=box, method=method) if(query.shape[0] == 3): - distances = mda.lib.distances.distance_array(query[None, :], - points, box=box) - else: - distances = mda.lib.distances.distance_array(query, + query = query.reshape((1, 3)) + + distances = mda.lib.distances.distance_array(query, points, box=box) - indices = np.where(distances < max_cutoff) + + if min_cutoff is None: + min_cutoff = 0. + + indices = np.where((distances < max_cutoff) & (distances > min_cutoff)) assert_equal(pairs[:, 1], indices[1]) From 7e25ea3ea1d0b88f9a1a296bc06f960f40553afb Mon Sep 17 00:00:00 2001 From: ayush Date: Thu, 12 Jul 2018 01:03:55 -0700 Subject: [PATCH 15/17] added few more tests --- .../MDAnalysisTests/lib/test_distances.py | 24 ++++++++++++------- 1 file changed, 16 insertions(+), 8 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index e9ff40f6777..acbf375c61b 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -54,13 +54,17 @@ def test_capped_distance_noresults(): assert_equal(len(pairs), 0) +npoints_1 = (1, 100) + boxes_1 = (np.array([1, 2, 3, 90, 90, 90], dtype=np.float32), # ortho np.array([1, 2, 3, 30, 45, 60], dtype=np.float32), # tri_box triclinic_vectors( # tri_vecs np.array([1, 2, 3, 90, 90, 45], dtype=np.float32)), np.array([[0.5, 0.9, 1.9], # tri_vecs_bad [2.0, 0.4, 0.1], - [0.0, 0.6, 0.5]], dtype=np.float32)) + [0.0, 0.6, 0.5]], dtype=np.float32), + None, # Non Periodic + ) query_1 = (np.array([0.1, 0.1, 0.1], dtype=np.float32), np.array([[0.1, 0.1, 0.1], @@ -70,15 +74,15 @@ def test_capped_distance_noresults(): min_cutoff_1 = (None, 0.1) - +@pytest.mark.parametrize('npoints', npoints_1) @pytest.mark.parametrize('box', boxes_1) @pytest.mark.parametrize('query', query_1) @pytest.mark.parametrize('method', method_1) @pytest.mark.parametrize('min_cutoff', min_cutoff_1) -def test_capped_distance_checkbrute(box, query, method, min_cutoff): +def test_capped_distance_checkbrute(npoints, box, query, method, min_cutoff): np.random.seed(90003) points = (np.random.uniform(low=0, high=1.0, - size=(100, 3))*(boxes_1[0][:3])).astype(np.float32) + size=(npoints, 3))*(boxes_1[0][:3])).astype(np.float32) max_cutoff = 0.3 # capped distance should be able to handle array of vectors # as well as single vectors. @@ -88,18 +92,22 @@ def test_capped_distance_checkbrute(box, query, method, min_cutoff): min_cutoff=min_cutoff, box=box, method=method) + if pairs.shape != (0, ): + found_pairs = pairs[:, 1] + else: + found_pairs = list() + if(query.shape[0] == 3): query = query.reshape((1, 3)) distances = mda.lib.distances.distance_array(query, - points, box=box) + points, box=box) if min_cutoff is None: min_cutoff = 0. - indices = np.where((distances < max_cutoff) & (distances > min_cutoff)) - - assert_equal(pairs[:, 1], indices[1]) + + assert_equal(np.sort(found_pairs, axis=0), np.sort(indices[1],axis=0)) # different boxlengths to shift a coordinate From edc5b2dcfb06350c6ea8bff36d7a3836aeabc864 Mon Sep 17 00:00:00 2001 From: ayush Date: Thu, 12 Jul 2018 01:58:13 -0700 Subject: [PATCH 16/17] more tests --- .../MDAnalysisTests/lib/test_distances.py | 44 ++++++++++++++++--- 1 file changed, 39 insertions(+), 5 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index acbf375c61b..f6e0dadf2cf 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -29,7 +29,7 @@ import MDAnalysis as mda -from MDAnalysis.lib.mdamath import triclinic_vectors +from MDAnalysis.lib.mdamath import triclinic_vectors, triclinic_box @pytest.mark.parametrize('coord_dtype', (np.float32, np.float64)) def test_transform_StoR_pass(coord_dtype): @@ -63,9 +63,10 @@ def test_capped_distance_noresults(): np.array([[0.5, 0.9, 1.9], # tri_vecs_bad [2.0, 0.4, 0.1], [0.0, 0.6, 0.5]], dtype=np.float32), - None, # Non Periodic + None, # Non Periodic ) + query_1 = (np.array([0.1, 0.1, 0.1], dtype=np.float32), np.array([[0.1, 0.1, 0.1], [0.2, 0.1, 0.1]], dtype=np.float32)) @@ -74,6 +75,7 @@ def test_capped_distance_noresults(): min_cutoff_1 = (None, 0.1) + @pytest.mark.parametrize('npoints', npoints_1) @pytest.mark.parametrize('box', boxes_1) @pytest.mark.parametrize('query', query_1) @@ -96,7 +98,7 @@ def test_capped_distance_checkbrute(npoints, box, query, method, min_cutoff): found_pairs = pairs[:, 1] else: found_pairs = list() - + if(query.shape[0] == 3): query = query.reshape((1, 3)) @@ -106,8 +108,40 @@ def test_capped_distance_checkbrute(npoints, box, query, method, min_cutoff): if min_cutoff is None: min_cutoff = 0. indices = np.where((distances < max_cutoff) & (distances > min_cutoff)) - - assert_equal(np.sort(found_pairs, axis=0), np.sort(indices[1],axis=0)) + + assert_equal(np.sort(found_pairs, axis=0), np.sort(indices[1], axis=0)) + + +npoints_2 = (1, 6000, 200000) + +cutoff_2 = (0.02, 0.2) + +methods = ( + '_bruteforce_capped', + '_bruteforce_capped', + '_pkdtree_capped', + '_pkdtree_capped', + '_pkdtree_capped', + '_bruteforce_capped', + ) + + +@pytest.mark.parametrize('ncm', zip(itertools.product(npoints_2, cutoff_2), methods)) +def test_method_selection(ncm): + np.random.seed(90003) + box = np.array([1, 1, 1, 90, 90, 90], dtype=np.float32) + points = (np.random.uniform(low=0, high=1.0, + size=(ncm[0][0], 3))*(box[:3])).astype(np.float32) + if box is not None: + boxtype = mda.lib.distances._box_check(box) + # Convert [A,B,C,alpha,beta,gamma] to [[A],[B],[C]] + if (boxtype == 'tri_box'): + box = triclinic_vectors(box) + if (boxtype == 'tri_vecs_bad'): + box = triclinic_vectors(triclinic_box(box[0], box[1], box[2])) + method = mda.lib.distances._determine_method(points, points, + ncm[0][1], box=box) + assert_equal(method.__name__, ncm[1]) # different boxlengths to shift a coordinate From 9d3a2fc15efd3f99e8f5bd62610a92b06b4fa1e0 Mon Sep 17 00:00:00 2001 From: Richard Gowers Date: Thu, 12 Jul 2018 10:06:35 -0500 Subject: [PATCH 17/17] Update test_distances.py --- .../MDAnalysisTests/lib/test_distances.py | 30 +++++++------------ 1 file changed, 11 insertions(+), 19 deletions(-) diff --git a/testsuite/MDAnalysisTests/lib/test_distances.py b/testsuite/MDAnalysisTests/lib/test_distances.py index f6e0dadf2cf..f61a30a2e37 100644 --- a/testsuite/MDAnalysisTests/lib/test_distances.py +++ b/testsuite/MDAnalysisTests/lib/test_distances.py @@ -112,26 +112,18 @@ def test_capped_distance_checkbrute(npoints, box, query, method, min_cutoff): assert_equal(np.sort(found_pairs, axis=0), np.sort(indices[1], axis=0)) -npoints_2 = (1, 6000, 200000) - -cutoff_2 = (0.02, 0.2) - -methods = ( - '_bruteforce_capped', - '_bruteforce_capped', - '_pkdtree_capped', - '_pkdtree_capped', - '_pkdtree_capped', - '_bruteforce_capped', - ) - - -@pytest.mark.parametrize('ncm', zip(itertools.product(npoints_2, cutoff_2), methods)) -def test_method_selection(ncm): +@pytest.mark.parametrize('npoints,cutoff,meth', + [(1, 0.02, '_bruteforce_capped'), + (1, 0.2, '_bruteforce_capped'), + (6000, 0.02, '_pkdtree_capped'), + (6000, 0.2, '_pkdtree_capped'), + (200000, 0.02, '_pkdtree_capped'), + (200000, 0.2, '_bruteforce_capped')]) +def test_method_selection(npoints, cutoff, meth): np.random.seed(90003) box = np.array([1, 1, 1, 90, 90, 90], dtype=np.float32) points = (np.random.uniform(low=0, high=1.0, - size=(ncm[0][0], 3))*(box[:3])).astype(np.float32) + size=(npoints, 3)) * (box[:3])).astype(np.float32) if box is not None: boxtype = mda.lib.distances._box_check(box) # Convert [A,B,C,alpha,beta,gamma] to [[A],[B],[C]] @@ -140,8 +132,8 @@ def test_method_selection(ncm): if (boxtype == 'tri_vecs_bad'): box = triclinic_vectors(triclinic_box(box[0], box[1], box[2])) method = mda.lib.distances._determine_method(points, points, - ncm[0][1], box=box) - assert_equal(method.__name__, ncm[1]) + cutoff, box=box) + assert_equal(method.__name__, meth) # different boxlengths to shift a coordinate