1#![allow(clippy::manual_div_ceil)]
2use std::fmt::Write as FmtWrite;
39
40use crate::Error as IoError;
41
42#[allow(dead_code)]
46#[derive(Debug, Clone, PartialEq)]
47pub struct FoamBoundaryPatch {
48 pub name: String,
50 pub patch_type: String,
52 pub uniform_value: Option<f64>,
54}
55
56pub fn foam_boundary_patch(
60 name: impl Into<String>,
61 patch_type: impl Into<String>,
62 uniform_value: Option<f64>,
63) -> FoamBoundaryPatch {
64 FoamBoundaryPatch {
65 name: name.into(),
66 patch_type: patch_type.into(),
67 uniform_value,
68 }
69}
70
71#[allow(dead_code)]
78#[derive(Debug, Clone)]
79pub struct FoamField {
80 pub name: String,
82 pub internal_field: Vec<f64>,
85 pub n_components: usize,
87 pub patches: Vec<FoamBoundaryPatch>,
89}
90
91impl FoamField {
92 pub fn new_scalar(name: impl Into<String>, values: Vec<f64>) -> Self {
94 Self {
95 name: name.into(),
96 internal_field: values,
97 n_components: 1,
98 patches: Vec::new(),
99 }
100 }
101
102 pub fn new_vector(name: impl Into<String>, values: Vec<f64>) -> Self {
104 Self {
105 name: name.into(),
106 internal_field: values,
107 n_components: 3,
108 patches: Vec::new(),
109 }
110 }
111
112 pub fn n_cells(&self) -> usize {
114 self.internal_field
115 .len()
116 .checked_div(self.n_components)
117 .unwrap_or(0)
118 }
119}
120
121fn extract_paren_block(s: &str) -> Option<&str> {
127 let start = s.find('(')?;
128 let end = s.rfind(')')?;
129 if end > start {
130 Some(&s[start + 1..end])
131 } else {
132 None
133 }
134}
135
136fn parse_scalar_list(list_text: &str) -> Result<Vec<f64>, IoError> {
138 let body = extract_paren_block(list_text)
139 .ok_or_else(|| IoError::Parse("missing ( ) in scalar list".into()))?;
140 body.split_whitespace()
141 .filter(|s| !s.is_empty())
142 .map(|tok| {
143 tok.parse::<f64>()
144 .map_err(|e| IoError::Parse(e.to_string()))
145 })
146 .collect()
147}
148
149fn parse_vector_list(list_text: &str) -> Result<Vec<f64>, IoError> {
151 let body = extract_paren_block(list_text)
152 .ok_or_else(|| IoError::Parse("missing outer ( ) in vector list".into()))?;
153 let mut values: Vec<f64> = Vec::new();
154 let mut remaining = body.trim();
155 while !remaining.is_empty() {
156 if let Some(open) = remaining.find('(') {
157 let close = remaining[open..].find(')').map(|i| open + i);
158 if let Some(close) = close {
159 let triple = &remaining[open + 1..close];
160 for tok in triple.split_whitespace() {
161 values.push(
162 tok.parse::<f64>()
163 .map_err(|e| IoError::Parse(e.to_string()))?,
164 );
165 }
166 remaining = &remaining[close + 1..];
167 } else {
168 break;
169 }
170 } else {
171 break;
172 }
173 }
174 Ok(values)
175}
176
177pub fn parse_foam_scalar(src: &str, name: impl Into<String>) -> Result<FoamField, IoError> {
186 let marker = "List<scalar>";
188 if let Some(pos) = src.find(marker) {
189 let after = &src[pos + marker.len()..];
190 let values = parse_scalar_list(after)?;
192 let mut field = FoamField::new_scalar(name, values);
193 field.patches = parse_boundary_patches(src);
195 return Ok(field);
196 }
197 let uniform_marker = "internalField uniform";
199 if let Some(pos) = src.find(uniform_marker) {
200 let after = &src[pos + uniform_marker.len()..];
201 let tok = after
202 .split_whitespace()
203 .next()
204 .ok_or_else(|| IoError::Parse("missing uniform value".into()))?;
205 let val: f64 = tok
206 .trim_end_matches(';')
207 .parse()
208 .map_err(|e: std::num::ParseFloatError| IoError::Parse(e.to_string()))?;
209 let n_cells = detect_n_cells(src);
210 let mut field = FoamField::new_scalar(name, vec![val; n_cells]);
211 field.patches = parse_boundary_patches(src);
212 return Ok(field);
213 }
214 Err(IoError::Parse(
215 "could not find internalField in foam scalar file".into(),
216 ))
217}
218
219pub fn parse_foam_vector(src: &str, name: impl Into<String>) -> Result<FoamField, IoError> {
226 let marker = "List<vector>";
227 if let Some(pos) = src.find(marker) {
228 let after = &src[pos + marker.len()..];
229 let values = parse_vector_list(after)?;
230 let mut field = FoamField::new_vector(name, values);
231 field.patches = parse_boundary_patches(src);
232 return Ok(field);
233 }
234 Err(IoError::Parse(
235 "could not find List<vector> in foam vector file".into(),
236 ))
237}
238
239fn detect_n_cells(src: &str) -> usize {
241 for line in src.lines() {
242 let trimmed = line.trim();
243 if let Some(rest) = trimmed.strip_prefix("// nCells:")
244 && let Ok(n) = rest.trim().parse::<usize>()
245 {
246 return n;
247 }
248 }
249 1
250}
251
252fn parse_boundary_patches(src: &str) -> Vec<FoamBoundaryPatch> {
256 let mut patches = Vec::new();
257 let bf_marker = "boundaryField";
259 let Some(bf_pos) = src.find(bf_marker) else {
260 return patches;
261 };
262 let after_bf = &src[bf_pos + bf_marker.len()..];
263 let Some(open) = after_bf.find('{') else {
265 return patches;
266 };
267 let inner = &after_bf[open + 1..];
268 let mut cursor = inner;
270 loop {
271 cursor = cursor.trim_start();
272 if cursor.is_empty() || cursor.starts_with('}') {
273 break;
274 }
275 let tok_end = cursor
277 .find(|c: char| c.is_whitespace() || c == '{')
278 .unwrap_or(cursor.len());
279 let patch_name = cursor[..tok_end].trim().to_string();
280 if patch_name.is_empty() || patch_name == "}" {
281 break;
282 }
283 cursor = &cursor[tok_end..];
285 let Some(blk_open) = cursor.find('{') else {
287 break;
288 };
289 let blk_start = blk_open + 1;
290 let Some(blk_close) = cursor[blk_start..].find('}') else {
291 break;
292 };
293 let block = &cursor[blk_start..blk_start + blk_close];
294 cursor = &cursor[blk_start + blk_close + 1..];
295 let patch_type = block
297 .lines()
298 .find_map(|l| {
299 let l = l.trim();
300 if l.starts_with("type") {
301 Some(
302 l.trim_start_matches("type")
303 .trim()
304 .trim_end_matches(';')
305 .trim()
306 .to_string(),
307 )
308 } else {
309 None
310 }
311 })
312 .unwrap_or_else(|| "unknown".to_string());
313 let uniform_value = block.lines().find_map(|l| {
315 let l = l.trim();
316 if l.starts_with("value") && l.contains("uniform") {
317 let parts: Vec<&str> = l.split_whitespace().collect();
318 if parts.len() >= 3 {
320 parts[2].trim_end_matches(';').parse::<f64>().ok()
321 } else {
322 None
323 }
324 } else {
325 None
326 }
327 });
328 patches.push(FoamBoundaryPatch {
329 name: patch_name,
330 patch_type,
331 uniform_value,
332 });
333 }
334 patches
335}
336
337pub fn write_foam_field(field: &FoamField) -> String {
339 let mut out = String::new();
340 let _ = writeln!(out, "// OpenFOAM field: {}", field.name);
341 let _ = writeln!(out, "dimensions [0 0 0 0 0 0 0];");
342 let _ = writeln!(out);
343 if field.n_components == 1 {
344 let _ = writeln!(out, "internalField nonuniform List<scalar>");
345 let _ = writeln!(out, "{}", field.n_cells());
346 let _ = writeln!(out, "(");
347 for &v in &field.internal_field {
348 let _ = writeln!(out, "{:.15e}", v);
349 }
350 let _ = writeln!(out, ");");
351 } else {
352 let _ = writeln!(out, "internalField nonuniform List<vector>");
353 let n = field.n_cells();
354 let _ = writeln!(out, "{}", n);
355 let _ = writeln!(out, "(");
356 for i in 0..n {
357 let base = i * 3;
358 let (vx, vy, vz) = (
359 field.internal_field[base],
360 field.internal_field[base + 1],
361 field.internal_field[base + 2],
362 );
363 let _ = writeln!(out, "({:.15e} {:.15e} {:.15e})", vx, vy, vz);
364 }
365 let _ = writeln!(out, ");");
366 }
367 let _ = writeln!(out);
368 let _ = writeln!(out, "boundaryField");
369 let _ = writeln!(out, "{{");
370 for p in &field.patches {
371 let _ = writeln!(out, " {}", p.name);
372 let _ = writeln!(out, " {{");
373 let _ = writeln!(out, " type {};", p.patch_type);
374 if let Some(val) = p.uniform_value {
375 let _ = writeln!(out, " value uniform {:.15e};", val);
376 }
377 let _ = writeln!(out, " }}");
378 }
379 let _ = writeln!(out, "}}");
380 out
381}
382
383#[allow(dead_code)]
390#[derive(Debug, Clone)]
391pub struct FoamMesh {
392 pub points: Vec<f64>,
394 pub faces: Vec<Vec<usize>>,
396 pub owner: Vec<usize>,
398 pub neighbour: Vec<usize>,
400}
401
402impl FoamMesh {
403 pub fn n_points(&self) -> usize {
405 self.points.len() / 3
406 }
407
408 pub fn n_faces(&self) -> usize {
410 self.faces.len()
411 }
412
413 pub fn n_internal_faces(&self) -> usize {
415 self.neighbour.len()
416 }
417}
418
419fn parse_foam_points(src: &str) -> Result<Vec<f64>, IoError> {
423 let body = extract_paren_block(src)
424 .ok_or_else(|| IoError::Parse("missing ( ) in points file".into()))?;
425 let mut pts: Vec<f64> = Vec::new();
426 for chunk in body.split('(').skip(1) {
427 let close = chunk
428 .find(')')
429 .ok_or_else(|| IoError::Parse("missing ) in point".into()))?;
430 let triple = &chunk[..close];
431 for tok in triple.split_whitespace() {
432 pts.push(
433 tok.parse::<f64>()
434 .map_err(|e| IoError::Parse(e.to_string()))?,
435 );
436 }
437 }
438 Ok(pts)
439}
440
441fn parse_foam_index_list(src: &str) -> Result<Vec<usize>, IoError> {
443 let body = extract_paren_block(src)
444 .ok_or_else(|| IoError::Parse("missing ( ) in index list".into()))?;
445 body.split_whitespace()
446 .filter(|s| !s.is_empty())
447 .map(|tok| {
448 tok.parse::<usize>()
449 .map_err(|e| IoError::Parse(e.to_string()))
450 })
451 .collect()
452}
453
454fn parse_foam_faces(src: &str) -> Result<Vec<Vec<usize>>, IoError> {
458 let body = extract_paren_block(src)
459 .ok_or_else(|| IoError::Parse("missing outer ( ) in faces file".into()))?;
460 let mut faces: Vec<Vec<usize>> = Vec::new();
461 let mut cursor = body.trim();
462 loop {
463 cursor = cursor.trim_start();
464 if cursor.is_empty() {
465 break;
466 }
467 let Some(open) = cursor.find('(') else { break };
469 let close = cursor[open..].find(')').map(|i| open + i);
470 let Some(close) = close else { break };
471 let face_str = &cursor[open + 1..close];
472 let pts: Vec<usize> = face_str
473 .split_whitespace()
474 .map(|t| {
475 t.parse::<usize>()
476 .map_err(|e| IoError::Parse(e.to_string()))
477 })
478 .collect::<Result<Vec<_>, _>>()?;
479 faces.push(pts);
480 cursor = &cursor[close + 1..];
481 }
482 Ok(faces)
483}
484
485pub fn read_foam_mesh(
493 points_src: &str,
494 faces_src: &str,
495 owner_src: &str,
496 neighbour_src: &str,
497) -> Result<FoamMesh, IoError> {
498 let points = parse_foam_points(points_src)?;
499 let faces = parse_foam_faces(faces_src)?;
500 let owner = parse_foam_index_list(owner_src)?;
501 let neighbour = parse_foam_index_list(neighbour_src)?;
502 Ok(FoamMesh {
503 points,
504 faces,
505 owner,
506 neighbour,
507 })
508}
509
510#[cfg(test)]
513mod tests {
514 use super::*;
515
516 const SCALAR_FILE: &str = r#"
517// OpenFOAM field
518dimensions [0 0 0 0 0 0 0];
519
520internalField nonuniform List<scalar>
5213
522(
5231.0
5242.5
5253.0
526);
527
528boundaryField
529{
530 inlet
531 {
532 type fixedValue;
533 value uniform 0;
534 }
535 outlet
536 {
537 type zeroGradient;
538 }
539}
540"#;
541
542 const VECTOR_FILE: &str = r#"
543// OpenFOAM vector field
544dimensions [0 1 -1 0 0 0 0];
545
546internalField nonuniform List<vector>
5472
548(
549(1.0 2.0 3.0)
550(4.0 5.0 6.0)
551);
552
553boundaryField
554{
555 wall
556 {
557 type noSlip;
558 }
559}
560"#;
561
562 const POINTS_FILE: &str = r#"
5633
564(
565(0.0 0.0 0.0)
566(1.0 0.0 0.0)
567(0.5 1.0 0.0)
568)
569"#;
570
571 const FACES_FILE: &str = r#"
5721
573(
5743(0 1 2)
575)
576"#;
577
578 const OWNER_FILE: &str = r#"
5791
580(
5810
582)
583"#;
584
585 const NEIGHBOUR_FILE: &str = r#"
5860
587(
588)
589"#;
590
591 #[test]
592 fn test_parse_scalar_n_cells() {
593 let f = parse_foam_scalar(SCALAR_FILE, "p").unwrap();
594 assert_eq!(f.n_cells(), 3);
595 }
596
597 #[test]
598 fn test_parse_scalar_values() {
599 let f = parse_foam_scalar(SCALAR_FILE, "p").unwrap();
600 assert!((f.internal_field[0] - 1.0).abs() < 1e-10);
601 assert!((f.internal_field[1] - 2.5).abs() < 1e-10);
602 assert!((f.internal_field[2] - 3.0).abs() < 1e-10);
603 }
604
605 #[test]
606 fn test_parse_scalar_name() {
607 let f = parse_foam_scalar(SCALAR_FILE, "pressure").unwrap();
608 assert_eq!(f.name, "pressure");
609 }
610
611 #[test]
612 fn test_parse_scalar_n_components() {
613 let f = parse_foam_scalar(SCALAR_FILE, "p").unwrap();
614 assert_eq!(f.n_components, 1);
615 }
616
617 #[test]
618 fn test_parse_scalar_patches() {
619 let f = parse_foam_scalar(SCALAR_FILE, "p").unwrap();
620 assert_eq!(f.patches.len(), 2);
621 assert_eq!(f.patches[0].name, "inlet");
622 assert_eq!(f.patches[0].patch_type, "fixedValue");
623 assert_eq!(f.patches[0].uniform_value, Some(0.0));
624 }
625
626 #[test]
627 fn test_parse_scalar_outlet_no_value() {
628 let f = parse_foam_scalar(SCALAR_FILE, "p").unwrap();
629 assert_eq!(f.patches[1].name, "outlet");
630 assert!(f.patches[1].uniform_value.is_none());
631 }
632
633 #[test]
634 fn test_parse_scalar_error_on_empty() {
635 assert!(parse_foam_scalar("", "p").is_err());
636 }
637
638 #[test]
639 fn test_parse_vector_n_cells() {
640 let f = parse_foam_vector(VECTOR_FILE, "U").unwrap();
641 assert_eq!(f.n_cells(), 2);
642 }
643
644 #[test]
645 fn test_parse_vector_n_components() {
646 let f = parse_foam_vector(VECTOR_FILE, "U").unwrap();
647 assert_eq!(f.n_components, 3);
648 }
649
650 #[test]
651 fn test_parse_vector_values() {
652 let f = parse_foam_vector(VECTOR_FILE, "U").unwrap();
653 assert!((f.internal_field[0] - 1.0).abs() < 1e-10);
654 assert!((f.internal_field[1] - 2.0).abs() < 1e-10);
655 assert!((f.internal_field[2] - 3.0).abs() < 1e-10);
656 assert!((f.internal_field[3] - 4.0).abs() < 1e-10);
657 }
658
659 #[test]
660 fn test_parse_vector_error_on_empty() {
661 assert!(parse_foam_vector("", "U").is_err());
662 }
663
664 #[test]
665 fn test_write_scalar_roundtrip() {
666 let original = FoamField {
667 name: "p".to_string(),
668 internal_field: vec![1.5, 2.5, 3.5],
669 n_components: 1,
670 patches: vec![foam_boundary_patch("inlet", "fixedValue", Some(0.0))],
671 };
672 let s = write_foam_field(&original);
673 let parsed = parse_foam_scalar(&s, "p").unwrap();
674 assert_eq!(parsed.n_cells(), 3);
675 assert!((parsed.internal_field[0] - 1.5).abs() < 1e-10);
676 assert!((parsed.internal_field[2] - 3.5).abs() < 1e-10);
677 }
678
679 #[test]
680 fn test_write_vector_roundtrip() {
681 let original = FoamField::new_vector("U", vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0]);
682 let s = write_foam_field(&original);
683 let parsed = parse_foam_vector(&s, "U").unwrap();
684 assert_eq!(parsed.n_cells(), 2);
685 for (i, &expected) in [1.0_f64, 2.0, 3.0, 4.0, 5.0, 6.0].iter().enumerate() {
686 assert!(
687 (parsed.internal_field[i] - expected).abs() < 1e-10,
688 "mismatch at {i}"
689 );
690 }
691 }
692
693 #[test]
694 fn test_write_field_contains_boundary_field() {
695 let f = FoamField::new_scalar("p", vec![0.0]);
696 let s = write_foam_field(&f);
697 assert!(s.contains("boundaryField"));
698 }
699
700 #[test]
701 fn test_foam_boundary_patch_constructor() {
702 let p = foam_boundary_patch("wall", "noSlip", None);
703 assert_eq!(p.name, "wall");
704 assert_eq!(p.patch_type, "noSlip");
705 assert!(p.uniform_value.is_none());
706 }
707
708 #[test]
709 fn test_foam_boundary_patch_with_value() {
710 let p = foam_boundary_patch("inlet", "fixedValue", Some(1.0));
711 assert_eq!(p.uniform_value, Some(1.0));
712 }
713
714 #[test]
715 fn test_read_foam_mesh_n_points() {
716 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
717 assert_eq!(m.n_points(), 3);
718 }
719
720 #[test]
721 fn test_read_foam_mesh_point_values() {
722 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
723 assert!((m.points[0] - 0.0).abs() < 1e-12);
724 assert!((m.points[3] - 1.0).abs() < 1e-12);
725 }
726
727 #[test]
728 fn test_read_foam_mesh_n_faces() {
729 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
730 assert_eq!(m.n_faces(), 1);
731 }
732
733 #[test]
734 fn test_read_foam_mesh_face_connectivity() {
735 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
736 assert_eq!(m.faces[0], vec![0, 1, 2]);
737 }
738
739 #[test]
740 fn test_read_foam_mesh_owner() {
741 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
742 assert_eq!(m.owner, vec![0]);
743 }
744
745 #[test]
746 fn test_read_foam_mesh_no_neighbours() {
747 let m = read_foam_mesh(POINTS_FILE, FACES_FILE, OWNER_FILE, NEIGHBOUR_FILE).unwrap();
748 assert_eq!(m.n_internal_faces(), 0);
749 }
750
751 #[test]
752 fn test_foam_field_new_scalar() {
753 let f = FoamField::new_scalar("T", vec![300.0, 310.0]);
754 assert_eq!(f.n_components, 1);
755 assert_eq!(f.n_cells(), 2);
756 }
757
758 #[test]
759 fn test_foam_field_new_vector() {
760 let f = FoamField::new_vector("U", vec![1.0, 0.0, 0.0, 2.0, 0.0, 0.0]);
761 assert_eq!(f.n_components, 3);
762 assert_eq!(f.n_cells(), 2);
763 }
764
765 #[test]
766 fn test_foam_field_debug() {
767 let f = FoamField::new_scalar("p", vec![]);
768 let s = format!("{f:?}");
769 assert!(s.contains("FoamField"));
770 }
771
772 #[test]
773 fn test_foam_mesh_debug() {
774 let m = FoamMesh {
775 points: vec![],
776 faces: vec![],
777 owner: vec![],
778 neighbour: vec![],
779 };
780 let s = format!("{m:?}");
781 assert!(s.contains("FoamMesh"));
782 }
783
784 #[test]
785 fn test_write_field_has_internal_field_keyword() {
786 let f = FoamField::new_scalar("rho", vec![1.0]);
787 let s = write_foam_field(&f);
788 assert!(s.contains("internalField"));
789 }
790
791 #[test]
792 fn test_parse_vector_patches() {
793 let f = parse_foam_vector(VECTOR_FILE, "U").unwrap();
794 assert_eq!(f.patches.len(), 1);
795 assert_eq!(f.patches[0].name, "wall");
796 }
797
798 #[test]
799 fn test_foam_boundary_patch_clone() {
800 let p = foam_boundary_patch("inlet", "fixedValue", Some(5.0));
801 let p2 = p.clone();
802 assert_eq!(p2.name, "inlet");
803 }
804
805 #[test]
806 fn test_foam_mesh_clone() {
807 let m = FoamMesh {
808 points: vec![0.0, 0.0, 0.0],
809 faces: vec![vec![0]],
810 owner: vec![0],
811 neighbour: vec![],
812 };
813 let m2 = m.clone();
814 assert_eq!(m2.n_points(), 1);
815 }
816}