use num_traits::AsPrimitive;
pub fn update_bvhnodes<Index>(bvhnodes: &mut [[Index; 3]], idx2vtx: &[Index], idx2morton: &[u32])
where
Index: num_traits::PrimInt + 'static + Copy,
usize: AsPrimitive<Index>,
{
assert_eq!(idx2vtx.len(), idx2morton.len());
assert!(!idx2morton.is_empty());
assert_eq!(bvhnodes.len(), idx2morton.len() * 2 - 1);
bvhnodes[0][0] = Index::max_value();
let num_branch = idx2morton.len() - 1; for i_branch in 0..num_branch {
let range = crate::mortons::range_of_binary_radix_tree_node(idx2morton, i_branch);
let isplit = crate::mortons::split_of_binary_radix_tree_node(idx2morton, range.0, range.1);
assert_ne!(isplit, usize::MAX);
if range.0 == isplit {
let i_bvhnode_a = num_branch + isplit; bvhnodes[i_branch][1] = i_bvhnode_a.as_();
bvhnodes[i_bvhnode_a][0] = i_branch.as_();
bvhnodes[i_bvhnode_a][1] = idx2vtx[isplit];
bvhnodes[i_bvhnode_a][2] = Index::max_value();
} else {
let i_bvhnode_a = isplit;
bvhnodes[i_branch][1] = i_bvhnode_a.as_();
bvhnodes[i_bvhnode_a][0] = i_branch.as_();
}
if range.1 == isplit + 1 {
let i_bvhnode_b = num_branch + isplit + 1;
bvhnodes[i_branch][2] = i_bvhnode_b.as_();
bvhnodes[i_bvhnode_b][0] = i_branch.as_();
bvhnodes[i_bvhnode_b][1] = idx2vtx[isplit + 1];
bvhnodes[i_bvhnode_b][2] = Index::max_value();
} else {
let i_bvhnode_b = isplit + 1;
bvhnodes[i_branch][2] = i_bvhnode_b.as_();
bvhnodes[i_bvhnode_b][0] = i_branch.as_();
}
}
}
#[test]
fn test_3d() {
let num_vtx = 100000;
let vtx2xyz: Vec<f32> = {
use rand::RngExt;
let mut rng = rand::rng();
(0..num_vtx * 3).map(|_| rng.random::<f32>()).collect()
};
let mut idx2vtx = vec![0usize; num_vtx];
let mut idx2morton = vec![0u32; num_vtx];
let mut vtx2morton = vec![0u32; num_vtx];
crate::mortons::sorted_morten_code3(
&mut idx2vtx,
&mut idx2morton,
&mut vtx2morton,
vtx2xyz.as_chunks::<3>().0,
&del_geo_core::mat4_col_major::from_identity(),
);
for idx in 0..num_vtx - 1 {
let jdx = idx + 1;
assert!(idx2morton[idx] <= idx2morton[jdx]);
}
crate::mortons::check_morton_code_range_split(&idx2morton);
let mut bvhnodes = vec![[0usize; 3]; num_vtx * 2 - 1];
update_bvhnodes(&mut bvhnodes, &idx2vtx, &idx2morton);
crate::bvhnodes::check_bvh_topology(&bvhnodes, num_vtx);
}
#[test]
fn test_2d() {
let num_vtx = 100000;
let vtx2xy: Vec<f32> = {
use rand::RngExt;
let mut rng = rand::rng();
(0..num_vtx * 2).map(|_| rng.random::<f32>()).collect()
};
let mut idx2vtx = vec![0usize; num_vtx];
let mut idx2morton = vec![0u32; num_vtx];
let mut vtx2morton = vec![0u32; num_vtx];
crate::mortons::sorted_morten_code2(
&mut idx2vtx,
&mut idx2morton,
&mut vtx2morton,
vtx2xy.as_chunks::<2>().0,
&del_geo_core::mat3_col_major::from_identity::<f32>(),
);
for idx in 0..num_vtx - 1 {
let jdx = idx + 1;
assert!(idx2morton[idx] <= idx2morton[jdx]);
}
crate::mortons::check_morton_code_range_split(&idx2morton);
let mut bvhnodes = vec![0usize; (num_vtx * 2 - 1) * 3];
update_bvhnodes(bvhnodes.as_chunks_mut::<3>().0, &idx2vtx, &idx2morton);
crate::bvhnodes::check_bvh_topology(bvhnodes.as_chunks::<3>().0, num_vtx);
}
pub fn from_vtx2xyz<Index, const NDIM: usize>(vtx2xyz: &[[f32; NDIM]]) -> Vec<[Index; 3]>
where
Index: num_traits::PrimInt + AsPrimitive<usize>,
usize: AsPrimitive<Index>,
{
let num_vtx = vtx2xyz.len();
let mut idx2vtx = vec![Index::one(); num_vtx];
let mut idx2morton = vec![0u32; num_vtx];
let mut vtx2morton = vec![0u32; num_vtx];
crate::mortons::update_sorted_morton_code::<_, NDIM>(
&mut idx2vtx,
&mut idx2morton,
&mut vtx2morton,
vtx2xyz,
);
let mut bvhnodes = vec![[Index::zero(); 3]; num_vtx * 2 - 1];
update_bvhnodes(&mut bvhnodes, &idx2vtx, &idx2morton);
bvhnodes
}
pub fn from_triangle_mesh<Index, const NDIM: usize>(
tri2vtx: &[[Index; 3]],
vtx2xy: &[[f32; NDIM]],
) -> Vec<[Index; 3]>
where
Index: num_traits::PrimInt + AsPrimitive<usize>,
usize: AsPrimitive<Index>,
{
let tri2cntr =
crate::elem2center::from_uniform_mesh_as_points::<Index, f32, 3, NDIM>(tri2vtx, vtx2xy);
from_vtx2xyz::<_, NDIM>(&tri2cntr)
}
pub fn update_for_triangle_mesh<Index>(
bvhnodes: &mut [[Index; 3]],
tri2vtx: &[[Index; 3]],
vtx2xy: &[f32],
num_dim: usize,
) where
Index: num_traits::PrimInt + AsPrimitive<usize>,
usize: AsPrimitive<Index>,
{
let num_tri = tri2vtx.len();
assert_eq!(bvhnodes.len(), num_tri * 2 - 1);
let mut idx2tri = vec![Index::one(); num_tri];
let mut idx2morton = vec![0u32; num_tri];
let mut tri2morton = vec![0u32; num_tri];
match num_dim {
2 => {
let tri2cntr = crate::elem2center::from_uniform_mesh_as_points::<Index, f32, 3, 2>(
tri2vtx,
vtx2xy.as_chunks::<2>().0,
);
crate::mortons::update_sorted_morton_code::<_, 2>(
&mut idx2tri,
&mut idx2morton,
&mut tri2morton,
&tri2cntr,
);
}
3 => {
let tri2cntr = crate::elem2center::from_uniform_mesh_as_points::<Index, f32, 3, 3>(
tri2vtx,
vtx2xy.as_chunks::<3>().0,
);
crate::mortons::update_sorted_morton_code::<_, 3>(
&mut idx2tri,
&mut idx2morton,
&mut tri2morton,
&tri2cntr,
);
}
_ => panic!("unsupported num_dim: {num_dim}"),
};
update_bvhnodes(bvhnodes, &idx2tri, &idx2morton);
}