From 7152218463a3a48fd18fde97b9f50eac67ac599e Mon Sep 17 00:00:00 2001 From: Dieter Weber Date: Fri, 22 Sep 2023 11:30:34 +0200 Subject: [PATCH 1/4] Faster implementation using ogrid When profiling the performance of diffraction pattern simulation, most time was spent in this function. This implementation runs more than three times faster since it reduces the calculation effort through broadcasting and re-uses coordinate and distance calculations. Instead of using diffpy, it performs the equivalent transformation to cartesian and calculation of the norm directly to allow these optimizations. --- diffsims/utils/sim_utils.py | 38 ++++++++++++++++++++++--------------- 1 file changed, 23 insertions(+), 15 deletions(-) diff --git a/diffsims/utils/sim_utils.py b/diffsims/utils/sim_utils.py index 09af2e80..5e5683fa 100644 --- a/diffsims/utils/sim_utils.py +++ b/diffsims/utils/sim_utils.py @@ -18,6 +18,7 @@ import collections import math +from functools import lru_cache import diffpy.structure import numpy as np @@ -447,21 +448,28 @@ def get_points_in_sphere(reciprocal_lattice, reciprocal_radius): Distance of reciprocal lattice points in sphere from the origin. """ a, b, c = reciprocal_lattice.a, reciprocal_lattice.b, reciprocal_lattice.c - h_max = np.floor(reciprocal_radius / a) - k_max = np.floor(reciprocal_radius / b) - l_max = np.floor(reciprocal_radius / c) - from itertools import product - - h_list = np.arange(-h_max, h_max + 1) # arange has a non-inclusive endpoint - k_list = np.arange(-k_max, k_max + 1) - l_list = np.arange(-l_max, l_max + 1) - potential_points = np.asarray(list(product(h_list, k_list, l_list))) - in_sphere = ( - np.abs(reciprocal_lattice.dist(potential_points, [0, 0, 0])) < reciprocal_radius - ) - spot_indices = potential_points[in_sphere] - cartesian_coordinates = reciprocal_lattice.cartesian(spot_indices) - spot_distances = reciprocal_lattice.dist(spot_indices, [0, 0, 0]) + h_max = int(np.floor(reciprocal_radius / a)) + k_max = int(np.floor(reciprocal_radius / b)) + l_max = int(np.floor(reciprocal_radius / c)) + + limit = reciprocal_radius**2 + h_list, k_list, l_list = np.ogrid[-h_max:h_max + 1, -k_max:k_max + 1, -l_max:l_max + 1] + full_shape = (2*h_max + 1, 2*k_max + 1, 2*l_max + 1) + + h_vec = h_list[..., np.newaxis] * reciprocal_lattice.base[0] + k_vec = k_list[..., np.newaxis] * reciprocal_lattice.base[1] + l_vec = l_list[..., np.newaxis] * reciprocal_lattice.base[2] + + vec = h_vec + k_vec + l_vec + norm2 = np.sum(vec**2, axis=-1) + select = norm2 < limit + h_select = np.broadcast_to(h_list, full_shape)[select] + k_select = np.broadcast_to(k_list, full_shape)[select] + l_select = np.broadcast_to(l_list, full_shape)[select] + spot_indices = np.stack((h_select, k_select, l_select), axis=-1) + + cartesian_coordinates = vec[select] + spot_distances = np.sqrt(norm2[select]) return spot_indices, cartesian_coordinates, spot_distances From 2505bc90a9f158b36407bf13f30b5fdf16b30e6a Mon Sep 17 00:00:00 2001 From: Dieter Weber Date: Fri, 22 Sep 2023 11:44:09 +0200 Subject: [PATCH 2/4] fix whitespace --- diffsims/utils/sim_utils.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/diffsims/utils/sim_utils.py b/diffsims/utils/sim_utils.py index 5e5683fa..166803cd 100644 --- a/diffsims/utils/sim_utils.py +++ b/diffsims/utils/sim_utils.py @@ -455,11 +455,11 @@ def get_points_in_sphere(reciprocal_lattice, reciprocal_radius): limit = reciprocal_radius**2 h_list, k_list, l_list = np.ogrid[-h_max:h_max + 1, -k_max:k_max + 1, -l_max:l_max + 1] full_shape = (2*h_max + 1, 2*k_max + 1, 2*l_max + 1) - + h_vec = h_list[..., np.newaxis] * reciprocal_lattice.base[0] k_vec = k_list[..., np.newaxis] * reciprocal_lattice.base[1] l_vec = l_list[..., np.newaxis] * reciprocal_lattice.base[2] - + vec = h_vec + k_vec + l_vec norm2 = np.sum(vec**2, axis=-1) select = norm2 < limit @@ -467,7 +467,7 @@ def get_points_in_sphere(reciprocal_lattice, reciprocal_radius): k_select = np.broadcast_to(k_list, full_shape)[select] l_select = np.broadcast_to(l_list, full_shape)[select] spot_indices = np.stack((h_select, k_select, l_select), axis=-1) - + cartesian_coordinates = vec[select] spot_distances = np.sqrt(norm2[select]) From c783f6f69901f4747d4ead52c5ad1235c548e3cd Mon Sep 17 00:00:00 2001 From: Dieter Weber Date: Fri, 22 Sep 2023 11:46:57 +0200 Subject: [PATCH 3/4] Add to credits, Zenodo --- .zenodo.json | 5 +++++ diffsims/release_info.py | 1 + 2 files changed, 6 insertions(+) diff --git a/.zenodo.json b/.zenodo.json index 6b19a039..8dcc119c 100644 --- a/.zenodo.json +++ b/.zenodo.json @@ -60,6 +60,11 @@ }, { "name":"Eirik Opheim" + }, + { + "name":"Dieter Weber", + "orcid": "0000-0001-6635-9567", + "affiliation": "Jülich Research Centre, Ernst Ruska Centre" } ] } diff --git a/diffsims/release_info.py b/diffsims/release_info.py index 4ac3c1c5..537e2beb 100644 --- a/diffsims/release_info.py +++ b/diffsims/release_info.py @@ -23,6 +23,7 @@ "Tina Bergh", "Tomas Ostasevicius", "Eirik Opheim", + "Dieter Weber", ] license = "GPLv3+" maintainer = "Duncan Johnstone, Phillip Crout, Håkon Wiik Ånes" From f072acc877634ca82477239befdb3526879843aa Mon Sep 17 00:00:00 2001 From: Dieter Weber Date: Fri, 22 Sep 2023 12:09:18 +0200 Subject: [PATCH 4/4] Remove unused import --- diffsims/utils/sim_utils.py | 1 - 1 file changed, 1 deletion(-) diff --git a/diffsims/utils/sim_utils.py b/diffsims/utils/sim_utils.py index 166803cd..fd2c414b 100644 --- a/diffsims/utils/sim_utils.py +++ b/diffsims/utils/sim_utils.py @@ -18,7 +18,6 @@ import collections import math -from functools import lru_cache import diffpy.structure import numpy as np