1use crate::offset::offset_surface;
28use crate::sweep_topology::parameter_line;
29use crate::topology::{
30 BrepSolid, CoedgeRecord, EdgeRecord, FaceRecord, LoopRecord, ShellRecord, VertexRecord,
31};
32use crate::image_curve::{affine_image_curve, image_curve_pair};
33use crate::{make_line, NurbsCurve, NurbsSurface, Vec3};
34
35fn principal_curvatures(surface: &NurbsSurface, u: f64, v: f64) -> Result<(f64, f64), String> {
38 surface
39 .principal_curvatures(u, v)
40 .map_err(|error| format!("thickenSheet: {error}"))
41}
42
43fn ensure_offsets_regular(surface: &NurbsSurface, distances: &[f64]) -> Result<(), String> {
48 let [u0, u1] = surface.domain_u()?;
49 let [v0, v1] = surface.domain_v()?;
50 const SAMPLES: usize = 33;
51 for i in 0..SAMPLES {
52 let u = u0 + (u1 - u0) * i as f64 / (SAMPLES - 1) as f64;
53 for j in 0..SAMPLES {
54 let v = v0 + (v1 - v0) * j as f64 / (SAMPLES - 1) as f64;
55 let (kappa_min, kappa_max) = principal_curvatures(surface, u, v)?;
56 for &distance in distances {
57 if distance == 0.0 {
58 continue;
59 }
60 for kappa in [kappa_min, kappa_max] {
61 let factor = 1.0 - distance * kappa;
62 if factor <= 1e-6 {
63 let radius = 1.0 / kappa.abs().max(1e-300);
64 return Err(format!(
65 "thickenSheet: offset by {distance:.6} self-intersects — the \
66 sheet's concave curvature radius {radius:.6} at (u={u:.4}, \
67 v={v:.4}) is not larger than the offset distance"
68 ));
69 }
70 }
71 }
72 }
73 }
74 Ok(())
75}
76
77fn offset_sheet(surface: &NurbsSurface, distance: f64) -> Result<NurbsSurface, String> {
81 if distance == 0.0 {
82 return Ok(surface.clone());
83 }
84 let carrier = FaceRecord {
85 id: 1,
86 surface: surface.clone(),
87 same_sense: true,
88 loops: vec![],
89 name: None,
90 };
91 offset_surface(&carrier, -distance, 0.0)
92}
93
94fn ruled_wall(bottom: &NurbsCurve, top: &NurbsCurve) -> Result<NurbsSurface, String> {
100 if bottom.degree != top.degree
101 || bottom.knots.len() != top.knots.len()
102 || bottom
103 .knots
104 .iter()
105 .zip(&top.knots)
106 .any(|(a, b)| (a - b).abs() > 1e-12)
107 || bottom
108 .control_points
109 .iter()
110 .zip(&top.control_points)
111 .any(|(a, b)| (a.w - b.w).abs() > 1e-9)
112 {
113 return Err("thickenSheet: internal error — offset sheet basis mismatch".into());
114 }
115 let rows = bottom
116 .control_points
117 .iter()
118 .zip(&top.control_points)
119 .map(|(b, t)| vec![*b, *t])
120 .collect();
121 NurbsSurface::new(
122 bottom.degree,
123 1,
124 bottom.knots.clone(),
125 vec![0.0, 0.0, 1.0, 1.0],
126 rows,
127 )
128}
129
130const GAUSS_X: [f64; 8] = [
131 -0.9602898564975363,
132 -0.7966664774136267,
133 -0.525532409916329,
134 -0.18343464249564978,
135 0.18343464249564978,
136 0.525532409916329,
137 0.7966664774136267,
138 0.9602898564975363,
139];
140const GAUSS_W: [f64; 8] = [
141 0.10122853629037669,
142 0.22238103445337445,
143 0.31370664587788727,
144 0.362683783378362,
145 0.362683783378362,
146 0.31370664587788727,
147 0.22238103445337445,
148 0.10122853629037669,
149];
150
151fn pcurve_signed_area(curve: &NurbsCurve) -> Result<f64, String> {
156 let [q0, q1] = curve.domain()?;
157 let mut breaks = vec![q0];
158 for &knot in &curve.knots {
159 if knot > q0 + 1e-12 && knot < q1 - 1e-12 && (knot - breaks[breaks.len() - 1]).abs() > 1e-12
160 {
161 breaks.push(knot);
162 }
163 }
164 breaks.push(q1);
165 let mut area = 0.0;
166 for pair in breaks.windows(2) {
167 let half = (pair[1] - pair[0]) * 0.5;
168 let middle = (pair[1] + pair[0]) * 0.5;
169 for index in 0..GAUSS_X.len() {
170 let derivatives = curve.derivatives(middle + half * GAUSS_X[index], 1)?;
171 let point = derivatives[0];
172 let tangent = derivatives[1];
173 area += GAUSS_W[index] * half * 0.5 * (point.x * tangent.y - point.y * tangent.x);
174 }
175 }
176 Ok(area)
177}
178
179fn planar_gap(first: Vec3, second: Vec3) -> f64 {
182 let du = first.x - second.x;
183 let dv = first.y - second.y;
184 (du * du + dv * dv).sqrt()
185}
186
187struct BoundaryImages {
191 bottom: NurbsCurve,
192 top: NurbsCurve,
193 t0: f64,
194 t1: f64,
195 dir: bool,
197}
198
199fn boundary_images(
216 base_affine: bool,
217 bottom: &NurbsSurface,
218 top: &NurbsSurface,
219 pcurve: &NurbsCurve,
220 eps_u: f64,
221 eps_v: f64,
222 fit_tolerance: f64,
223) -> Result<BoundaryImages, String> {
224 if base_affine {
225 let [q0, q1] = pcurve.domain()?;
226 return Ok(BoundaryImages {
227 bottom: affine_image_curve(bottom, pcurve)?,
228 top: affine_image_curve(top, pcurve)?,
229 t0: q0,
230 t1: q1,
231 dir: true,
232 });
233 }
234 if pcurve.degree == 1 && pcurve.control_points.len() == 2 {
235 let first = pcurve.control_points[0];
236 let second = pcurve.control_points[1];
237 if (first.w - second.w).abs() <= 1e-12 {
238 let (ua, va) = (first.x / first.w, first.y / first.w);
239 let (ub, vb) = (second.x / second.w, second.y / second.w);
240 if (ua - ub).abs() <= eps_u && (va - vb).abs() > eps_v {
241 let u_constant = (ua + ub) * 0.5;
242 return Ok(BoundaryImages {
243 bottom: bottom.iso_curve_u(u_constant)?,
244 top: top.iso_curve_u(u_constant)?,
245 t0: va.min(vb),
246 t1: va.max(vb),
247 dir: vb > va,
248 });
249 }
250 if (va - vb).abs() <= eps_v && (ua - ub).abs() > eps_u {
251 let v_constant = (va + vb) * 0.5;
252 return Ok(BoundaryImages {
253 bottom: bottom.iso_curve_v(v_constant)?,
254 top: top.iso_curve_v(v_constant)?,
255 t0: ua.min(ub),
256 t1: ua.max(ub),
257 dir: ub > ua,
258 });
259 }
260 }
261 }
262 let (bottom_image, top_image) =
265 image_curve_pair(bottom, top, pcurve, fit_tolerance, "thickenSheet")?;
266 let forward = bottom_image.t0 <= bottom_image.t1;
267 Ok(BoundaryImages {
268 bottom: bottom_image.curve,
269 top: top_image.curve,
270 t0: bottom_image.t0.min(bottom_image.t1),
271 t1: bottom_image.t0.max(bottom_image.t1),
272 dir: forward,
273 })
274}
275
276pub fn thicken_trimmed_sheet(
313 surface: &NurbsSurface,
314 loops: &[Vec<NurbsCurve>],
315 thickness: f64,
316 symmetric: bool,
317) -> Result<BrepSolid, String> {
318 if !thickness.is_finite() || thickness.abs() <= 1e-12 {
319 return Err("thickenSheet: thickness must be a nonzero finite value".into());
320 }
321 if loops.is_empty() || loops.iter().any(|loop_curves| loop_curves.is_empty()) {
322 return Err("thickenSheet: at least one non-empty pcurve loop is required".into());
323 }
324 let (closed_u, closed_v) = surface.closed_directions()?;
325 if closed_u || closed_v {
326 return Err(
327 "thickenSheet: closed sheets are not supported (split the patch at its seam first)"
328 .into(),
329 );
330 }
331 let (distance_bottom, distance_top) = if symmetric {
332 (-thickness.abs() * 0.5, thickness.abs() * 0.5)
333 } else if thickness > 0.0 {
334 (0.0, thickness)
335 } else {
336 (thickness, 0.0)
337 };
338 ensure_offsets_regular(surface, &[distance_bottom, distance_top])?;
339
340 let bottom = offset_sheet(surface, distance_bottom)?;
341 let top = offset_sheet(surface, distance_top)?;
342 let [u0, u1] = surface.domain_u()?;
343 let [v0, v1] = surface.domain_v()?;
344 let uv_tolerance = 1e-7 * (u1 - u0).max(v1 - v0);
345 let minimum_area = 1e-10 * (u1 - u0) * (v1 - v0);
346 let eps_u = 1e-9 * (u1 - u0);
347 let eps_v = 1e-9 * (v1 - v0);
348 let base_affine = surface.is_affine()?;
349 let sheet_points = bottom
356 .control_points
357 .iter()
358 .flatten()
359 .map(|control| control.point())
360 .collect::<Result<Vec<_>, String>>()?;
361 let fit_tolerance =
362 crate::KernelTolerances::for_scale(crate::model_scale(sheet_points), 1e-7).intersection_fit;
363
364 let mut vertices: Vec<VertexRecord> = Vec::new();
365 let mut edges: Vec<EdgeRecord> = Vec::new();
366 let mut faces: Vec<FaceRecord> = Vec::new();
367 let mut top_cap_loops: Vec<LoopRecord> = Vec::new();
368 let mut bottom_cap_loops: Vec<LoopRecord> = Vec::new();
369 let mut bottom_junction_points: Vec<Vec3> = Vec::new();
370 let mut next_id = 1u64;
371
372 for (loop_index, loop_curves) in loops.iter().enumerate() {
373 let count = loop_curves.len();
374
375 let mut starts = Vec::with_capacity(count);
377 let mut ends = Vec::with_capacity(count);
378 for curve in loop_curves {
379 let [q0, q1] = curve.domain()?;
380 starts.push(curve.evaluate(q0)?);
381 ends.push(curve.evaluate(q1)?);
382 }
383 for index in 0..count {
384 let next_index = (index + 1) % count;
385 let gap = planar_gap(ends[index], starts[next_index]);
386 if gap > uv_tolerance {
387 return Err(format!(
388 "thickenSheet: loop {loop_index} is open — pcurve {index} ends at \
389 (u={:.6}, v={:.6}) but pcurve {next_index} starts at (u={:.6}, v={:.6}) \
390 (parameter-space gap {gap:.3e})",
391 ends[index].x, ends[index].y, starts[next_index].x, starts[next_index].y
392 ));
393 }
394 }
395 if count == 1 {
396 let [q0, q1] = loop_curves[0].domain()?;
397 let middle = loop_curves[0].evaluate((q0 + q1) * 0.5)?;
398 if planar_gap(middle, starts[0]) <= uv_tolerance {
399 return Err(format!(
400 "thickenSheet: loop {loop_index} is degenerate (zero parameter-space extent)"
401 ));
402 }
403 } else {
404 for index in 0..count {
405 if planar_gap(ends[index], starts[index]) <= uv_tolerance {
406 return Err(format!(
407 "thickenSheet: pcurve {index} of loop {loop_index} closes on itself \
408 inside a multi-curve loop (pinched loop)"
409 ));
410 }
411 }
412 }
413 let mut area = 0.0;
414 for curve in loop_curves {
415 area += pcurve_signed_area(curve)?;
416 }
417 if loop_index == 0 {
418 if area <= minimum_area {
419 return Err(format!(
420 "thickenSheet: outer loop must run counter-clockwise in (u, v) \
421 (signed area {area:.3e})"
422 ));
423 }
424 } else if area >= -minimum_area {
425 return Err(format!(
426 "thickenSheet: hole loop {loop_index} must run clockwise in (u, v) \
427 (signed area {area:.3e})"
428 ));
429 }
430
431 let mut bottom_vertex_ids = Vec::with_capacity(count);
433 let mut top_vertex_ids = Vec::with_capacity(count);
434 let mut bottom_points = Vec::with_capacity(count);
435 let mut top_points = Vec::with_capacity(count);
436 for start in &starts {
437 let bottom_point = bottom.evaluate(start.x, start.y)?;
438 let top_point = top.evaluate(start.x, start.y)?;
439 vertices.push(VertexRecord {
440 id: next_id,
441 point: bottom_point,
442 });
443 bottom_vertex_ids.push(next_id);
444 next_id += 1;
445 vertices.push(VertexRecord {
446 id: next_id,
447 point: top_point,
448 });
449 top_vertex_ids.push(next_id);
450 next_id += 1;
451 bottom_points.push(bottom_point);
452 top_points.push(top_point);
453 bottom_junction_points.push(bottom_point);
454 }
455
456 let mut images = Vec::with_capacity(count);
458 for curve in loop_curves {
459 images.push(boundary_images(
460 base_affine,
461 &bottom,
462 &top,
463 curve,
464 eps_u,
465 eps_v,
466 fit_tolerance,
467 )?);
468 }
469 let mut bottom_edge_ids = Vec::with_capacity(count);
470 let mut top_edge_ids = Vec::with_capacity(count);
471 for (index, image) in images.iter().enumerate() {
472 let next_index = (index + 1) % count;
473 let (start_j, end_j) = if image.dir {
474 (index, next_index)
475 } else {
476 (next_index, index)
477 };
478 edges.push(EdgeRecord {
479 id: next_id,
480 curve: image.bottom.clone(),
481 t0: image.t0,
482 t1: image.t1,
483 start_vertex_id: bottom_vertex_ids[start_j],
484 end_vertex_id: bottom_vertex_ids[end_j],
485 degenerate: false,
486 name: None,
487 });
488 bottom_edge_ids.push(next_id);
489 next_id += 1;
490 edges.push(EdgeRecord {
491 id: next_id,
492 curve: image.top.clone(),
493 t0: image.t0,
494 t1: image.t1,
495 start_vertex_id: top_vertex_ids[start_j],
496 end_vertex_id: top_vertex_ids[end_j],
497 degenerate: false,
498 name: None,
499 });
500 top_edge_ids.push(next_id);
501 next_id += 1;
502 }
503 let mut vertical_edge_ids = Vec::with_capacity(count);
504 for junction in 0..count {
505 edges.push(EdgeRecord {
506 id: next_id,
507 curve: make_line(bottom_points[junction], top_points[junction])?,
508 t0: 0.0,
509 t1: 1.0,
510 start_vertex_id: bottom_vertex_ids[junction],
511 end_vertex_id: top_vertex_ids[junction],
512 degenerate: false,
513 name: None,
514 });
515 vertical_edge_ids.push(next_id);
516 next_id += 1;
517 }
518
519 for (index, image) in images.iter().enumerate() {
529 let next_index = (index + 1) % count;
530 let wall = ruled_wall(&image.bottom, &image.top)?;
531 let (s_start, s_end) = if image.dir {
532 (image.t0, image.t1)
533 } else {
534 (image.t1, image.t0)
535 };
536 let mut coedges = Vec::with_capacity(4);
537 for (edge_id, forward, pcurve) in [
538 (
539 bottom_edge_ids[index],
540 image.dir,
541 parameter_line(s_start, 0.0, s_end, 0.0)?,
542 ),
543 (
544 vertical_edge_ids[next_index],
545 true,
546 parameter_line(s_end, 0.0, s_end, 1.0)?,
547 ),
548 (
549 top_edge_ids[index],
550 !image.dir,
551 parameter_line(s_end, 1.0, s_start, 1.0)?,
552 ),
553 (
554 vertical_edge_ids[index],
555 false,
556 parameter_line(s_start, 1.0, s_start, 0.0)?,
557 ),
558 ] {
559 coedges.push(CoedgeRecord {
560 id: next_id,
561 edge_id,
562 forward,
563 pcurve,
564 });
565 next_id += 1;
566 }
567 let loop_id = next_id;
568 next_id += 1;
569 faces.push(FaceRecord {
570 id: next_id,
571 surface: wall,
572 same_sense: image.dir,
573 loops: vec![LoopRecord {
574 id: loop_id,
575 coedges,
576 }],
577 name: None,
578 });
579 next_id += 1;
580 }
581
582 let mut top_coedges = Vec::with_capacity(count);
584 for (index, image) in images.iter().enumerate() {
585 top_coedges.push(CoedgeRecord {
586 id: next_id,
587 edge_id: top_edge_ids[index],
588 forward: image.dir,
589 pcurve: loop_curves[index].clone(),
590 });
591 next_id += 1;
592 }
593 top_cap_loops.push(LoopRecord {
594 id: next_id,
595 coedges: top_coedges,
596 });
597 next_id += 1;
598 let mut bottom_coedges = Vec::with_capacity(count);
599 for index in (0..count).rev() {
600 bottom_coedges.push(CoedgeRecord {
601 id: next_id,
602 edge_id: bottom_edge_ids[index],
603 forward: !images[index].dir,
604 pcurve: loop_curves[index].reversed()?,
605 });
606 next_id += 1;
607 }
608 bottom_cap_loops.push(LoopRecord {
609 id: next_id,
610 coedges: bottom_coedges,
611 });
612 next_id += 1;
613 }
614
615 let scale = crate::model_scale(bottom_junction_points.iter().copied());
623 for first in 0..bottom_junction_points.len() {
624 for second in first + 1..bottom_junction_points.len() {
625 if bottom_junction_points[first]
626 .sub(bottom_junction_points[second])
627 .length()
628 <= 1e-7 * scale
629 {
630 return Err(
631 "thickenSheet: sheet boundary is degenerate (coincident junction vertices)"
632 .into(),
633 );
634 }
635 }
636 }
637
638 faces.push(FaceRecord {
641 id: next_id,
642 surface: top,
643 same_sense: true,
644 loops: top_cap_loops,
645 name: None,
646 });
647 next_id += 1;
648 faces.push(FaceRecord {
649 id: next_id,
650 surface: bottom,
651 same_sense: false,
652 loops: bottom_cap_loops,
653 name: None,
654 });
655 next_id += 1;
656
657 let shell_id = next_id;
658 let solid = BrepSolid {
659 id: next_id + 1,
660 vertices,
661 edges,
662 shells: vec![ShellRecord {
663 id: shell_id,
664 faces,
665 }],
666 genus: loops.len() as i64 - 1,
667 };
668 let issues = solid.validate();
669 if !issues.is_empty() {
670 return Err(format!(
671 "thickenSheet: assembled solid failed validation: {issues:?}"
672 ));
673 }
674 let volume = crate::solid_signed_volume(&solid)?;
675 if volume <= 0.0 {
676 return Err(format!(
677 "thickenSheet: internal orientation error (signed volume {volume})"
678 ));
679 }
680 Ok(solid)
681}
682
683pub fn thicken_face_sheet(
694 surface: &NurbsSurface,
695 thickness: f64,
696 symmetric: bool,
697) -> Result<BrepSolid, String> {
698 if !thickness.is_finite() || thickness.abs() <= 1e-12 {
699 return Err("thickenSheet: thickness must be a nonzero finite value".into());
700 }
701 let [u0, u1] = surface.domain_u()?;
702 let [v0, v1] = surface.domain_v()?;
703 let rectangle = vec![
704 parameter_line(u0, v0, u1, v0)?,
705 parameter_line(u1, v0, u1, v1)?,
706 parameter_line(u1, v1, u0, v1)?,
707 parameter_line(u0, v1, u0, v0)?,
708 ];
709 thicken_trimmed_sheet(surface, &[rectangle], thickness, symmetric)
710}
711
712