honeycomb-kernels 0.11.0

Meshing kernel implementations using combinatorial maps
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117

//! Step 5 implementation
//
//! Use the information computed at step 4 and insert all new edges into the map.

use honeycomb_core::cmap::{CMap2, DartIdType};
use honeycomb_core::geometry::CoordsFloat;
use honeycomb_core::stm::atomically_with_err;

use crate::cell_insertion::insert_vertices_on_edge;
use crate::grisubal::model::{Boundary, MapEdge};
use crate::utils::VertexAnchor;

#[allow(clippy::cast_possible_truncation)]
pub(crate) fn insert_edges_in_map<T: CoordsFloat>(cmap: &mut CMap2<T>, edges: &[MapEdge<T>]) {
    let dart_slices = build_workload(cmap, edges);

    for (
        i,
        (
            MapEdge {
                start,
                intermediates,
                end,
            },
            dslice,
        ),
    ) in edges.iter().zip(dart_slices.iter()).enumerate()
    {
        let &[d_new, b2_d_new] = &dslice[0..2] else {
            unreachable!()
        };
        build_base_edge(cmap, *start, *end, [d_new, b2_d_new]);

        if !intermediates.is_empty() {
            // create the topology components
            let edge_id = cmap.edge_id(d_new);
            let new_darts = &dslice[2..];
            atomically_with_err(|t| {
                insert_vertices_on_edge(
                    cmap,
                    t,
                    edge_id,
                    new_darts,
                    &vec![T::from(0.5).unwrap(); intermediates.len()],
                )
            })
            .unwrap();
            // replace placeholder vertices
            let mut dart_id = cmap.beta::<1>(edge_id as DartIdType);
            for v in intermediates {
                let vid = cmap.vertex_id(dart_id);
                let _ = cmap.write_vertex(vid, *v);
                if cmap.contains_attribute::<VertexAnchor>() {
                    let _ = cmap.write_attribute::<VertexAnchor>(vid, VertexAnchor::Node(i as u32));
                }
                dart_id = cmap.beta::<1>(dart_id);
            }
        }

        mark_boundary(cmap, *start, *end);
    }
}

#[allow(clippy::cast_possible_truncation)]
fn build_workload<T: CoordsFloat>(
    cmap: &mut CMap2<T>,
    edges: &[MapEdge<T>],
) -> Vec<Vec<DartIdType>> {
    // allocate all darts needed
    let n_tot: usize = edges.iter().map(|e| 2 + 2 * e.intermediates.len()).sum();
    let tmp = cmap.allocate_used_darts(n_tot) as usize;

    // the prefix sum gives an offset that corresponds to the starting index of each slice, minus
    // the location of the allocated dart block (given by `tmp`)
    // we can deduce theend of the slice from start index (prefix sum) + number of darts
    edges
        .iter()
        .map(|e| 2 + 2 * e.intermediates.len())
        .scan(0, |state, n_d| {
            *state += n_d;
            Some((n_d, *state - n_d)) // we want an offset, not the actual sum
        })
        .map(|(n_d, start)| {
            ((tmp + start) as DartIdType..(tmp + start + n_d) as DartIdType).collect::<Vec<_>>()
        })
        .collect()
}

fn build_base_edge<T: CoordsFloat>(
    cmap: &CMap2<T>,
    start: DartIdType,
    end: DartIdType,
    [d_new, b2_d_new]: [DartIdType; 2],
) {
    // remove deprecated connectivities & save what data is necessary
    let b1_start_old = cmap.beta::<1>(start);
    let b0_end_old = cmap.beta::<0>(end);
    cmap.unlink::<1>(start).unwrap();
    cmap.unlink::<1>(b0_end_old).unwrap();

    cmap.link::<2>(d_new, b2_d_new).unwrap();

    // rebuild - this is the final construct if there are no intermediates
    cmap.link::<1>(start, d_new).unwrap();
    cmap.link::<1>(b2_d_new, b1_start_old).unwrap();
    cmap.link::<1>(d_new, end).unwrap();
    cmap.link::<1>(b0_end_old, b2_d_new).unwrap();
}

fn mark_boundary<T: CoordsFloat>(cmap: &CMap2<T>, start: DartIdType, end: DartIdType) {
    let mut d_boundary = cmap.beta::<1>(start);
    while d_boundary != end {
        cmap.write_attribute::<Boundary>(d_boundary, Boundary::Left);
        cmap.write_attribute::<Boundary>(cmap.beta::<2>(d_boundary), Boundary::Right);
        d_boundary = cmap.beta::<1>(d_boundary);
    }
}