Created
June 6, 2024 22:13
-
-
Save akaszynski/da3295444738fcfcd8f2d56084b40f69 to your computer and use it in GitHub Desktop.
Demonstrate mapping face centers from an external source to the cells of an unstructured grid.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
""" | |
Demonstrate mapping face centers from an external source to the cells of an unstructured grid. | |
Assumptions: | |
- Points from the external source are close to the centers of the faces of the | |
extracted surface | |
""" | |
from pykdtree.kdtree import KDTree | |
import pyvista as pv | |
import numpy as np | |
def make_ugrid(density=30): | |
"""Construct an unstructured grid from a structured grid.""" | |
s = 1 / (density - 1) | |
ugrid = pv.ImageData( | |
dimensions=(density, density, density), spacing=(s, s, s) | |
).cast_to_unstructured_grid() | |
ugrid.points = ugrid.points.astype(np.float64) | |
# uniform grids generate voxels and we need hexahedrons cells | |
ugrid.celltypes[:] = pv.CellType.HEXAHEDRON | |
ccon = ugrid.cell_connectivity.reshape(-1, 8) | |
ccon[:, 7], ccon[:, 6] = ccon[:, 6], ccon[:, 7].copy() | |
ccon[:, 3], ccon[:, 2] = ccon[:, 2], ccon[:, 3].copy() | |
return ugrid | |
ugrid = make_ugrid() | |
surf = ugrid.extract_surface() | |
# apply a fiticious boundry condition for each "face" in our extracted surface | |
bc = np.zeros((surf.n_cells, 3)) | |
mask = surf.cell_centers().points[:, 2] == 0 | |
bc[mask, 2] = 1 | |
surf.cell_data["bc"] = bc | |
# for each cell in the grid | |
face_to_cell = surf["vtkOriginalCellIds"] | |
bc_grid = np.zeros((ugrid.n_cells, 3)) | |
bc_grid[face_to_cell] = bc | |
# Generate simulated face centers from an external source | |
ccent = surf.cell_centers().points | |
ccent += np.random.random(ccent.shape) * 1e-4 - 5e-4 | |
np.random.shuffle(ccent) # inplace | |
external_data = np.random.random(ccent.shape) | |
pl = pv.Plotter() | |
pl.add_title("Show external points overlapping the unstructured grid") | |
pl.add_points(ccent, point_size=10, color="r", render_points_as_spheres=True) | |
pl.add_mesh(ugrid, color="w") | |
pl.show() | |
# get nearest point of our surf for each point from the external data | |
kd_tree = KDTree(surf.cell_centers().points) | |
dist, ind = kd_tree.query(ccent, k=1) | |
# You now have a map between the external points to our surface. If you had | |
# data associated with each point of the external surface points (for example, | |
# `external_data`), you would then: | |
# | |
surf.cell_data["external_data"] = external_data[ind] | |
# if you needed this data mapped back to `ugrid` cells, then: | |
# | |
# create a zeroed array (because each cell won't have a surface face) | |
ext_ugrid_data = np.zeros((ugrid.n_cells, 3)) | |
ext_ugrid_data[face_to_cell] = surf.cell_data["external_data"] | |
ugrid["external_data"] = ext_ugrid_data | |
# now, each cell of ugrid contains the external data, though note that some | |
# information will be lost since each cell will contain multiple faces | |
############################################################################### | |
from pykdtree.kdtree import KDTree | |
import pyvista as pv | |
import numpy as np | |
ugrid = make_ugrid() | |
# goal: identify correspondence between external surface data and a cells within ugrid | |
# Generate simulated face centers from an external source | |
surf = ugrid.extract_surface() | |
ccent = surf.cell_centers().points | |
ccent += np.random.random(ccent.shape) * 1e-3 - 5e-4 | |
np.random.shuffle(ccent) # inplace | |
kd_tree = KDTree(ugrid.cell_centers().points) | |
dist, ind = kd_tree.query(ccent, k=1) | |
# For fun: each ind now corresponds to a cell center in ugrid | |
ugrid.extract_cells(ind).explode(2).plot() |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment