honeycomb_benches/
cut_edges.rs

1use std::time::Instant;
2
3use rayon::prelude::*;
4
5use honeycomb::{
6    core::{
7        cmap::SewError,
8        stm::{Transaction, TransactionControl, TransactionResult},
9    },
10    kernels::remeshing::{cut_inner_edge, cut_outer_edge},
11    prelude::{CMap2, CMapBuilder, CoordsFloat, DartIdType, EdgeIdType},
12};
13
14use crate::{cli::CutEdgesArgs, utils::hash_file};
15
16// const MAX_RETRY: u8 = 10;
17
18pub fn bench_cut_edges<T: CoordsFloat>(args: CutEdgesArgs) -> CMap2<T> {
19    let input_map = args.input.to_str().unwrap();
20    let target_len = T::from(args.target_length).unwrap();
21
22    let n_threads = std::thread::available_parallelism()
23        .map(|n| n.get())
24        .unwrap_or(1);
25
26    // load map from file
27    let mut instant = Instant::now();
28    let input_hash = hash_file(input_map).expect("E: could not compute input hash"); // file id for posterity
29
30    let mut map: CMap2<T> = if input_map.ends_with(".cmap") {
31        CMapBuilder::<2, _>::from_cmap_file(input_map)
32            .build()
33            .unwrap()
34    } else if input_map.ends_with(".vtk") {
35        CMapBuilder::<2, _>::from_vtk_file(input_map)
36            .build()
37            .unwrap()
38    } else {
39        panic!(
40            "E: Unknown file format; only .cmap or .vtk files are supported for map initialization"
41        );
42    };
43    #[cfg(debug_assertions)] // check input
44    {
45        use honeycomb::prelude::OrbitPolicy;
46        assert!(
47            map.iter_faces()
48                .all(|f| { map.orbit(OrbitPolicy::Face, f as DartIdType).count() == 3 }),
49            "Input mesh isn't a triangle mesh"
50        );
51    }
52    println!("| cut-edges benchmark");
53    println!("|-> input      : {input_map} (hash: {input_hash:#0x})",);
54    println!(
55        "|-> backend    : {:?} with {n_threads} thread(s)",
56        args.backend
57    );
58    println!("|-> target size: {target_len:?}");
59    println!("|-> init time  : {}ms", instant.elapsed().as_millis());
60
61    println!(
62        " Step | n_edge_total | n_edge_to_process | t_compute_batch(s) | t_process_batch(s) | n_transac_retry"
63    );
64
65    let mut step = 0;
66    print!(" {step:>4} "); // Step
67
68    // compute first batch
69    instant = Instant::now();
70    let mut edges: Vec<EdgeIdType> = map.iter_edges().collect();
71    print!("| {:>12} ", edges.len()); // n_edge_total
72    edges.retain(|&e| {
73        let (vid1, vid2) = (
74            map.vertex_id(e as DartIdType),
75            map.vertex_id(map.beta::<1>(e as DartIdType)),
76        );
77        match (map.force_read_vertex(vid1), map.force_read_vertex(vid2)) {
78            (Some(v1), Some(v2)) => (v2 - v1).norm() > target_len,
79            (_, _) => false,
80        }
81    });
82    let n_e = edges.len();
83    print!("| {n_e:>17} "); // n_edge_to_process
84    let mut nd = map.add_free_darts(6 * n_e); // 2 for edge split + 2*2 for new edges in neighbor tets
85    let mut darts: Vec<DartIdType> = (nd..nd + 6 * n_e as DartIdType).collect();
86    print!("| {:>18.6e} ", instant.elapsed().as_secs_f64()); // t_compute_batch
87
88    // while there are edges to cut
89    while !edges.is_empty() {
90        // process batch
91        instant = Instant::now();
92        let n_retry = match args.backend {
93            crate::cli::Backend::RayonIter => dispatch_rayon(&map, &mut edges, &darts),
94            crate::cli::Backend::RayonChunks => {
95                dispatch_rayon_chunks(&map, &mut edges, &darts, n_threads)
96            }
97            crate::cli::Backend::StdThreads => {
98                dispatch_std_threads(&map, &mut edges, &darts, n_threads)
99            }
100        };
101        print!("| {:>18.6e} ", instant.elapsed().as_secs_f64()); // t_process_batch
102        println!("| {n_retry:>15}",); // n_transac_retry
103
104        (1..map.n_darts() as DartIdType).for_each(|d| {
105            if map.is_free(d) && !map.is_unused(d) {
106                map.remove_free_dart(d);
107            }
108        });
109
110        // compute the new batch
111        step += 1;
112        print!(" {step:>4} "); // Step
113        instant = Instant::now();
114        edges.extend(map.iter_edges());
115        print!("| {:>12} ", edges.len()); // n_edge_total
116        edges.retain(|&e| {
117            let (vid1, vid2) = (
118                map.vertex_id(e as DartIdType),
119                map.vertex_id(map.beta::<1>(e as DartIdType)),
120            );
121            match (map.force_read_vertex(vid1), map.force_read_vertex(vid2)) {
122                (Some(v1), Some(v2)) => (v2 - v1).norm() > target_len,
123                (_, _) => false,
124            }
125        });
126        let n_e = edges.len();
127        print!("| {n_e:>17} "); // n_edge_to_process
128        nd = map.add_free_darts(6 * n_e);
129        darts.par_drain(..); // is there a better way?
130        darts.extend(nd..nd + 6 * n_e as DartIdType);
131        if n_e != 0 {
132            print!("| {:>18.6e} ", instant.elapsed().as_secs_f64()); // t_compute_batch
133        } else {
134            print!("| {:>18.6e} ", instant.elapsed().as_secs_f64()); // t_compute_batch
135            print!("| {:>18.6e} ", 0.0); // t_process_batch
136            println!("| {:>15}", 0); // n_transac_retry
137        }
138    }
139
140    map
141}
142
143#[inline]
144fn dispatch_rayon<T: CoordsFloat>(
145    map: &CMap2<T>,
146    edges: &mut Vec<EdgeIdType>,
147    darts: &[DartIdType],
148) -> u32 {
149    let units: Vec<(u32, [u32; 6])> = edges
150        .drain(..)
151        .zip(darts.chunks(6))
152        .map(|(e, sl)| (e, sl.try_into().unwrap()))
153        .collect();
154    units
155        .into_par_iter()
156        .map(|(e, new_darts)| {
157            let mut n_retry = 0;
158            if map.is_i_free::<2>(e as DartIdType) {
159                while !process_outer_edge(map, &mut n_retry, e, new_darts).is_validated() {}
160            } else {
161                while !process_inner_edge(map, &mut n_retry, e, new_darts).is_validated() {}
162            }
163            n_retry as u32
164        }) // par_map
165        .sum()
166}
167
168#[inline]
169fn dispatch_rayon_chunks<T: CoordsFloat>(
170    map: &CMap2<T>,
171    edges: &mut Vec<EdgeIdType>,
172    darts: &[DartIdType],
173    n_threads: usize,
174) -> u32 {
175    let units: Vec<(u32, [u32; 6])> = edges
176        .drain(..)
177        .zip(darts.chunks(6))
178        .map(|(e, sl)| (e, sl.try_into().unwrap()))
179        .collect();
180    units
181        .par_chunks(1 + units.len() / n_threads)
182        .map(|c| {
183            let mut n = 0;
184            c.iter().for_each(|&(e, new_darts)| {
185                let mut n_retry = 0;
186                if map.is_i_free::<2>(e as DartIdType) {
187                    while !process_outer_edge(map, &mut n_retry, e, new_darts).is_validated() {}
188                } else {
189                    while !process_inner_edge(map, &mut n_retry, e, new_darts).is_validated() {}
190                }
191                n += n_retry as u32;
192            });
193            n
194        }) // par_for_each
195        .sum()
196}
197
198#[inline]
199fn dispatch_std_threads<T: CoordsFloat>(
200    map: &CMap2<T>,
201    edges: &mut Vec<EdgeIdType>,
202    darts: &[DartIdType],
203    n_threads: usize,
204) -> u32 {
205    let units: Vec<(u32, [u32; 6])> = edges
206        .drain(..)
207        .zip(darts.chunks(6))
208        .map(|(e, sl)| (e, sl.try_into().unwrap()))
209        .collect();
210    std::thread::scope(|s| {
211        let mut handles = Vec::new();
212        for wl in units.chunks(1 + units.len() / n_threads) {
213            handles.push(s.spawn(|| {
214                let mut n = 0;
215                wl.iter().for_each(|&(e, new_darts)| {
216                    let mut n_retry = 0;
217                    if map.is_i_free::<2>(e as DartIdType) {
218                        while !process_outer_edge(map, &mut n_retry, e, new_darts).is_validated() {}
219                    } else {
220                        while !process_inner_edge(map, &mut n_retry, e, new_darts).is_validated() {}
221                    }
222                    n += n_retry as u32;
223                });
224                n
225            })); // s.spawn
226        } // for wl in workloads
227        handles.into_iter().map(|h| h.join().unwrap()).sum()
228    }) // std::thread::scope
229}
230
231#[inline]
232fn process_outer_edge<T: CoordsFloat>(
233    map: &CMap2<T>,
234    n_retry: &mut u8,
235    e: EdgeIdType,
236    [nd1, nd2, nd3, _, _, _]: [DartIdType; 6],
237) -> TransactionResult<(), SewError> {
238    Transaction::with_control_and_err(
239        |_| {
240            *n_retry += 1;
241            TransactionControl::Retry
242        },
243        |trans| cut_outer_edge(trans, map, e, [nd1, nd2, nd3]),
244    ) // Transaction::with_control
245}
246
247#[inline]
248fn process_inner_edge<T: CoordsFloat>(
249    map: &CMap2<T>,
250    n_retry: &mut u8,
251    e: EdgeIdType,
252    nds: [DartIdType; 6],
253) -> TransactionResult<(), SewError> {
254    Transaction::with_control_and_err(
255        |_| {
256            *n_retry += 1;
257            TransactionControl::Retry
258        },
259        |trans| cut_inner_edge(trans, map, e, nds),
260    ) // Transaction::with_control
261}
honeycomb_benches/cut_edges.rs

honeycomb_benches/
cut_edges.rs
