Skip to main content

runmat_meshing_size/
field.rs

1use serde::{Deserialize, Serialize};
2
3pub const MODULE_PURPOSE: &str = "composable sizing queries for every meshing stage";
4
5#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
6pub struct SizingSample {
7    pub position_m: [f64; 3],
8    pub target_size_m: f64,
9    #[serde(default)]
10    pub reason: Option<String>,
11}
12
13#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
14pub struct AnisotropicSizingSample {
15    pub position_m: [f64; 3],
16    pub target_sizes_m: [f64; 3],
17    pub directions: [[f64; 3]; 3],
18    #[serde(default)]
19    pub reason: Option<String>,
20}
21
22impl AnisotropicSizingSample {
23    pub fn is_valid_metric(&self) -> bool {
24        self.position_m.iter().all(|value| value.is_finite())
25            && self
26                .target_sizes_m
27                .iter()
28                .all(|value| value.is_finite() && *value > 0.0)
29            && directions_are_finite_orthonormal(self.directions)
30    }
31}
32
33#[derive(Debug, Clone, Copy, PartialEq)]
34pub struct SegmentSizingQuery {
35    pub start_m: [f64; 3],
36    pub end_m: [f64; 3],
37}
38
39#[derive(Debug, Clone, Copy, PartialEq)]
40pub struct PointSizingQuery {
41    pub position_m: [f64; 3],
42}
43
44#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
45#[serde(rename_all = "snake_case")]
46pub enum SizingQuerySource {
47    Unset,
48    Global,
49    LocalSample,
50    AnisotropicMetric,
51}
52
53#[derive(Debug, Clone, Copy, PartialEq, Serialize, Deserialize)]
54pub struct SizingQueryResult {
55    pub target_size_m: Option<f64>,
56    pub source: SizingQuerySource,
57    pub contributing_sample_count: usize,
58}
59
60impl SizingQueryResult {
61    pub fn unset() -> Self {
62        Self {
63            target_size_m: None,
64            source: SizingQuerySource::Unset,
65            contributing_sample_count: 0,
66        }
67    }
68}
69
70pub trait SizingFieldService {
71    fn query_point_size(&self, query: PointSizingQuery) -> SizingQueryResult;
72    fn query_segment_size(&self, query: SegmentSizingQuery) -> SizingQueryResult;
73}
74
75#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
76pub struct SizingSampleRejection {
77    pub position_m: [f64; 3],
78    pub target_size_m: f64,
79    pub status: String,
80    #[serde(default)]
81    pub reason: Option<String>,
82    #[serde(default)]
83    pub detail: Option<String>,
84}
85
86#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
87pub struct SizingSampleApplication {
88    pub position_m: [f64; 3],
89    pub target_size_m: f64,
90    pub inserted_breakpoint_count: usize,
91    #[serde(default)]
92    pub reason: Option<String>,
93    #[serde(default)]
94    pub detail: Option<String>,
95}
96
97#[derive(Debug, Clone, PartialEq, Serialize, Deserialize, Default)]
98pub struct MeshSizingField {
99    #[serde(default)]
100    pub global_target_size_m: Option<f64>,
101    #[serde(default)]
102    pub min_size_m: Option<f64>,
103    #[serde(default)]
104    pub max_size_m: Option<f64>,
105    #[serde(default)]
106    pub growth_rate: Option<f64>,
107    #[serde(default)]
108    pub samples: Vec<SizingSample>,
109    #[serde(default)]
110    pub anisotropic_samples: Vec<AnisotropicSizingSample>,
111    #[serde(default)]
112    pub applied_samples: Vec<SizingSampleApplication>,
113    #[serde(default)]
114    pub rejected_samples: Vec<SizingSampleRejection>,
115}
116
117impl MeshSizingField {
118    pub fn target_size_for_segment(&self, query: SegmentSizingQuery) -> Option<f64> {
119        self.query_segment_size(query).target_size_m
120    }
121
122    pub fn target_size_at_point(&self, query: PointSizingQuery) -> Option<f64> {
123        self.query_point_size(query).target_size_m
124    }
125
126    fn initial_query_result(&self) -> SizingQueryResult {
127        self.global_target_size_m
128            .and_then(|target_size_m| self.clamped_target_size_m(target_size_m))
129            .map(|target_size_m| SizingQueryResult {
130                target_size_m: Some(target_size_m),
131                source: SizingQuerySource::Global,
132                contributing_sample_count: 0,
133            })
134            .unwrap_or_else(SizingQueryResult::unset)
135    }
136
137    fn merge_query_target(
138        &self,
139        mut result: SizingQueryResult,
140        target_size_m: f64,
141        source: SizingQuerySource,
142    ) -> SizingQueryResult {
143        let Some(target_size_m) = self.clamped_target_size_m(target_size_m) else {
144            return result;
145        };
146        if match result.target_size_m {
147            Some(current) => target_size_m < current,
148            None => true,
149        } {
150            result.target_size_m = Some(target_size_m);
151            result.source = source;
152        }
153        result.contributing_sample_count += 1;
154        result
155    }
156
157    pub fn clamped_target_size_m(&self, target_size_m: f64) -> Option<f64> {
158        if !target_size_m.is_finite() || target_size_m <= 0.0 {
159            return None;
160        }
161        let mut target_size_m = target_size_m;
162        if let (Some(global_target_size_m), Some(growth_rate)) = (
163            self.global_target_size_m
164                .filter(|value| value.is_finite() && *value > 0.0),
165            self.growth_rate
166                .filter(|value| value.is_finite() && *value >= 1.0),
167        ) {
168            target_size_m = target_size_m.max(global_target_size_m / growth_rate);
169        }
170        if let Some(min_size_m) = self
171            .min_size_m
172            .filter(|value| value.is_finite() && *value > 0.0)
173        {
174            target_size_m = target_size_m.max(min_size_m);
175        }
176        if let Some(max_size_m) = self
177            .max_size_m
178            .filter(|value| value.is_finite() && *value > 0.0)
179        {
180            target_size_m = target_size_m.min(max_size_m);
181        }
182        (target_size_m.is_finite() && target_size_m > 0.0).then_some(target_size_m)
183    }
184}
185
186impl SizingFieldService for MeshSizingField {
187    fn query_point_size(&self, query: PointSizingQuery) -> SizingQueryResult {
188        if !query.position_m.iter().all(|value| value.is_finite()) {
189            return SizingQueryResult::unset();
190        }
191        let mut result = self.initial_query_result();
192        for sample in &self.samples {
193            if point_matches(sample.position_m, query.position_m) {
194                result = self.merge_query_target(
195                    result,
196                    sample.target_size_m,
197                    SizingQuerySource::LocalSample,
198                );
199            }
200        }
201        for sample in &self.anisotropic_samples {
202            if !sample.is_valid_metric() || !point_matches(sample.position_m, query.position_m) {
203                continue;
204            }
205            let sample_target_size_m = sample
206                .target_sizes_m
207                .iter()
208                .copied()
209                .fold(f64::INFINITY, f64::min);
210            result = self.merge_query_target(
211                result,
212                sample_target_size_m,
213                SizingQuerySource::AnisotropicMetric,
214            );
215        }
216        result
217    }
218
219    fn query_segment_size(&self, query: SegmentSizingQuery) -> SizingQueryResult {
220        let mut result = self.initial_query_result();
221        for sample in &self.samples {
222            if point_influences_segment(sample.position_m, sample.target_size_m, query, self) {
223                result = self.merge_query_target(
224                    result,
225                    sample.target_size_m,
226                    SizingQuerySource::LocalSample,
227                );
228            }
229        }
230        for sample in &self.anisotropic_samples {
231            let sample_target_size_m = sample
232                .target_sizes_m
233                .iter()
234                .copied()
235                .fold(f64::INFINITY, f64::min);
236            if !sample.is_valid_metric()
237                || !point_influences_segment(sample.position_m, sample_target_size_m, query, self)
238            {
239                continue;
240            }
241            result = self.merge_query_target(
242                result,
243                sample_target_size_m,
244                SizingQuerySource::AnisotropicMetric,
245            );
246        }
247        result
248    }
249}
250
251fn point_influences_segment(
252    point: [f64; 3],
253    target_size_m: f64,
254    query: SegmentSizingQuery,
255    sizing: &MeshSizingField,
256) -> bool {
257    if !point.iter().all(|value| value.is_finite())
258        || !query.start_m.iter().all(|value| value.is_finite())
259        || !query.end_m.iter().all(|value| value.is_finite())
260    {
261        return false;
262    }
263    let segment = sub(query.end_m, query.start_m);
264    let segment_length_squared = dot(segment, segment);
265    if segment_length_squared <= f64::EPSILON {
266        return distance(point, query.start_m) <= 1.0e-12;
267    }
268    let relative = sub(point, query.start_m);
269    let parameter = dot(relative, segment) / segment_length_squared;
270    if !(-1.0e-12..=1.0 + 1.0e-12).contains(&parameter) {
271        return false;
272    }
273    let closest = [
274        query.start_m[0] + segment[0] * parameter.clamp(0.0, 1.0),
275        query.start_m[1] + segment[1] * parameter.clamp(0.0, 1.0),
276        query.start_m[2] + segment[2] * parameter.clamp(0.0, 1.0),
277    ];
278    let tolerance = segment_length_squared.sqrt().max(1.0) * 1.0e-9;
279    let influence_radius = sizing
280        .clamped_target_size_m(target_size_m)
281        .unwrap_or(0.0)
282        .max(tolerance)
283        .max(1.0e-12);
284    distance(point, closest) <= influence_radius
285}
286
287fn point_matches(left: [f64; 3], right: [f64; 3]) -> bool {
288    if !left.iter().all(|value| value.is_finite()) || !right.iter().all(|value| value.is_finite()) {
289        return false;
290    }
291    let scale = left
292        .iter()
293        .chain(right.iter())
294        .map(|value| value.abs())
295        .fold(1.0_f64, f64::max);
296    distance(left, right) <= scale * 1.0e-9
297}
298
299fn sub(left: [f64; 3], right: [f64; 3]) -> [f64; 3] {
300    [left[0] - right[0], left[1] - right[1], left[2] - right[2]]
301}
302
303fn distance(left: [f64; 3], right: [f64; 3]) -> f64 {
304    dot(sub(left, right), sub(left, right)).sqrt()
305}
306
307fn directions_are_finite_orthonormal(directions: [[f64; 3]; 3]) -> bool {
308    directions
309        .iter()
310        .all(|direction| direction.iter().all(|value| value.is_finite()))
311        && directions.iter().all(|direction| {
312            let norm_squared = dot(*direction, *direction);
313            (norm_squared - 1.0).abs() <= 1.0e-9
314        })
315        && dot(directions[0], directions[1]).abs() <= 1.0e-9
316        && dot(directions[0], directions[2]).abs() <= 1.0e-9
317        && dot(directions[1], directions[2]).abs() <= 1.0e-9
318}
319
320fn dot(left: [f64; 3], right: [f64; 3]) -> f64 {
321    left[0] * right[0] + left[1] * right[1] + left[2] * right[2]
322}
323
324#[cfg(test)]
325mod tests {
326    use super::*;
327
328    #[test]
329    fn anisotropic_sizing_sample_validates_metric_contract() {
330        let valid = AnisotropicSizingSample {
331            position_m: [0.1, 0.2, 0.3],
332            target_sizes_m: [0.01, 0.02, 0.04],
333            directions: [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]],
334            reason: Some("boundary_layer".to_string()),
335        };
336        assert!(valid.is_valid_metric());
337
338        let mut invalid_size = valid.clone();
339        invalid_size.target_sizes_m[1] = 0.0;
340        assert!(!invalid_size.is_valid_metric());
341
342        let mut non_orthogonal = valid;
343        non_orthogonal.directions[1] = [1.0, 0.0, 0.0];
344        assert!(!non_orthogonal.is_valid_metric());
345    }
346
347    #[test]
348    fn mesh_sizing_field_deserializes_without_anisotropic_samples() {
349        let sizing: MeshSizingField =
350            serde_json::from_str(r#"{"global_target_size_m":0.1}"#).expect("sizing field");
351        assert_eq!(sizing.global_target_size_m, Some(0.1));
352        assert!(sizing.anisotropic_samples.is_empty());
353    }
354
355    #[test]
356    fn segment_sizing_query_uses_samples_on_segment() {
357        let sizing = MeshSizingField {
358            global_target_size_m: Some(1.0),
359            samples: vec![
360                SizingSample {
361                    position_m: [0.5, 0.0, 0.0],
362                    target_size_m: 0.2,
363                    reason: Some("feature_edge".to_string()),
364                },
365                SizingSample {
366                    position_m: [0.5, 0.2, 0.0],
367                    target_size_m: 0.05,
368                    reason: Some("off_edge".to_string()),
369                },
370            ],
371            ..MeshSizingField::default()
372        };
373
374        let target_size_m = sizing
375            .target_size_for_segment(SegmentSizingQuery {
376                start_m: [0.0, 0.0, 0.0],
377                end_m: [1.0, 0.0, 0.0],
378            })
379            .expect("segment target");
380
381        assert_eq!(target_size_m, 0.2);
382        let query = sizing.query_segment_size(SegmentSizingQuery {
383            start_m: [0.0, 0.0, 0.0],
384            end_m: [1.0, 0.0, 0.0],
385        });
386        assert_eq!(query.target_size_m, Some(0.2));
387        assert_eq!(query.source, SizingQuerySource::LocalSample);
388        assert_eq!(query.contributing_sample_count, 1);
389    }
390
391    #[test]
392    fn segment_sizing_query_uses_samples_within_target_radius() {
393        let sizing = MeshSizingField {
394            global_target_size_m: Some(1.0),
395            samples: vec![SizingSample {
396                position_m: [0.5, 0.25, 0.0],
397                target_size_m: 0.3,
398                reason: Some("boundary_patch".to_string()),
399            }],
400            ..MeshSizingField::default()
401        };
402
403        let query = sizing.query_segment_size(SegmentSizingQuery {
404            start_m: [0.0, 0.0, 0.0],
405            end_m: [1.0, 0.0, 0.0],
406        });
407
408        assert_eq!(query.target_size_m, Some(0.3));
409        assert_eq!(query.source, SizingQuerySource::LocalSample);
410        assert_eq!(query.contributing_sample_count, 1);
411    }
412
413    #[test]
414    fn segment_sizing_query_clamps_local_samples() {
415        let sizing = MeshSizingField {
416            global_target_size_m: Some(1.0),
417            min_size_m: Some(0.25),
418            max_size_m: Some(0.75),
419            samples: vec![SizingSample {
420                position_m: [0.5, 0.0, 0.0],
421                target_size_m: 0.1,
422                reason: Some("feature_edge".to_string()),
423            }],
424            ..MeshSizingField::default()
425        };
426
427        let target_size_m = sizing
428            .target_size_for_segment(SegmentSizingQuery {
429                start_m: [0.0, 0.0, 0.0],
430                end_m: [1.0, 0.0, 0.0],
431            })
432            .expect("segment target");
433
434        assert_eq!(target_size_m, 0.25);
435    }
436
437    #[test]
438    fn point_sizing_query_reports_anisotropic_metric_source() {
439        let sizing = MeshSizingField {
440            global_target_size_m: Some(1.0),
441            samples: vec![SizingSample {
442                position_m: [0.25, 0.0, 0.0],
443                target_size_m: 0.4,
444                reason: Some("feature".to_string()),
445            }],
446            anisotropic_samples: vec![AnisotropicSizingSample {
447                position_m: [0.25, 0.0, 0.0],
448                target_sizes_m: [0.2, 0.3, 0.5],
449                directions: [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]],
450                reason: Some("directional".to_string()),
451            }],
452            ..MeshSizingField::default()
453        };
454
455        let query = sizing.query_point_size(PointSizingQuery {
456            position_m: [0.25, 0.0, 0.0],
457        });
458
459        assert_eq!(query.target_size_m, Some(0.2));
460        assert_eq!(query.source, SizingQuerySource::AnisotropicMetric);
461        assert_eq!(query.contributing_sample_count, 2);
462        assert_eq!(
463            sizing.target_size_at_point(PointSizingQuery {
464                position_m: [0.25, 0.0, 0.0],
465            }),
466            Some(0.2)
467        );
468    }
469
470    #[test]
471    fn invalid_point_sizing_query_is_unset() {
472        let sizing = MeshSizingField {
473            global_target_size_m: Some(1.0),
474            ..MeshSizingField::default()
475        };
476
477        let query = sizing.query_point_size(PointSizingQuery {
478            position_m: [f64::NAN, 0.0, 0.0],
479        });
480
481        assert_eq!(query, SizingQueryResult::unset());
482    }
483}