From 2dd08ababcde3d901b3747f75a6ddcd7e4247ac5 Mon Sep 17 00:00:00 2001 From: Timotej Bernat Date: Thu, 25 Apr 2024 20:01:55 -0600 Subject: [PATCH] Dismantled coordops, distributed functionality amongst .integral and .bravais --- polymerist/maths/lattices/coordops.py | 66 --------------------------- 1 file changed, 66 deletions(-) delete mode 100644 polymerist/maths/lattices/coordops.py diff --git a/polymerist/maths/lattices/coordops.py b/polymerist/maths/lattices/coordops.py deleted file mode 100644 index 4487d11..0000000 --- a/polymerist/maths/lattices/coordops.py +++ /dev/null @@ -1,66 +0,0 @@ -'''Distance geometry methods for arrays of coordinates''' - -from typing import TypeVar -from ...genutils.typetools.numpytypes import Shape, N, M, DType -from ...genutils.typetools.categorical import Numeric -TT = TypeVar('TT', bound=Numeric) - -import numpy as np -from .coordinates import BoundingBox -from .integer import CubicIntegerLattice - - -def nearest_int_coord_along_normal(point : np.ndarray[Shape[N], Numeric], normal : np.ndarray[Shape[N], Numeric]) -> np.ndarray[Shape[N], int]: - ''' - Takes an N-dimensional control point and an N-dimensional normal vector - Returns the integer-valued point nearest to the control point which lies in - the normal direction relative to the control point - ''' - # Check that inputs have vector-like shapes and compatible dimensions - assert(point.ndim == 1) - assert(normal.ndim == 1) - assert(point.size == normal.size) - - min_int_bound_point = np.floor(point) - max_int_bound_point = np.ceil( point) - if np.isclose(min_int_bound_point, max_int_bound_point).all(): # edge case: if already on an integer-valued point, the fllor and ceiling will be identical - int_point = point - else: - int_bbox = BoundingBox(np.vstack([min_int_bound_point, max_int_bound_point])) # generate bounding box from extremal positions - dots = np.inner((int_bbox.vertices - point), normal) # take dot product between the normal and the direction vectors from the control point to the integer bounding points - i = dots.argmax() # position of integer point in most similar direction to normal - furthest_point, similarity = int_bbox.vertices[i], dots[i] - - if similarity <= 0: - raise ValueError(f'Could not locate valid integral point in normal direction (closest found was {furthest_point} with dot product {similarity})') - else: - int_point = furthest_point - - return int_point.astype(int) - -def identify_bravais_points_within_bbox(lattice_vectors : np.ndarray[Shape[N, N], TT], bbox : BoundingBox) -> tuple[np.ndarray[Shape[M, N], TT], CubicIntegerLattice]: - '''Locate all lattice points generated by a set of Bravais lattice vectors in N-dimensions which fall within a given bounding box - Returns the coordinate vector of the internal lattice points and a CubicIntegerLattice containing the corresponding lattice vector multiplicities''' - # 1) transform the bounding box into the inverse domain, where Bravais lattice points "look like" integer points - inv_bbox_coords = bbox.linear_transformation(np.linalg.inv(lattice_vectors.T), as_coords=True) - directors = (inv_bbox_coords.points - inv_bbox_coords.centroid) # directional vectors at each vertex which are away from transformed box center - directors /= np.linalg.norm(directors, axis=1)[:, None] # normalize direction vectors - - # 2) find the nearest integer-valued points to transformed bounding box points in directions "away" from the bounding box - integral_bounding_points = np.vstack([ - nearest_int_coord_along_normal(point, director) # ...(ensures no bounded points end up outside of final box) - for point, director in zip(inv_bbox_coords.points, directors) - ]) - - # 3) produce collection of integer point which contains AT LEAST the integral points bounded by the transformed nearest-integrer bounding points - trial_bbox = BoundingBox(integral_bounding_points) # TODO: box is pretty inefficient for highly-sheared unit vectors: need simple method to test for integral points within a bounding volume (or simplex) - trial_int_latt = CubicIntegerLattice(trial_bbox.dims + 1) # need K+1 unit-spaced points to space an interval K units long (e.g. [0-2] -> (0, 1, 2)) - trial_int_latt.points += trial_bbox.minimum # offset lattice from default at origin - - # 4) transforms the trial points back to the normal domain and performa final enclosure check - trial_latt_points = trial_int_latt.linear_transformation(lattice_vectors, as_coords=True) - is_inside = bbox.surrounds(trial_latt_points).all(axis=1) # mask of which points have ALL coordinates inside the bounding box - inside_latt_points = trial_latt_points.points[is_inside] - trial_int_latt.points = trial_int_latt.points[ is_inside] # screen out integral points which fall outside the bounding box when mapped to lattice points - - return inside_latt_points, trial_int_latt \ No newline at end of file