embedded3 PERF-3 P1-1: the moving path's host work, digit-identical — impose_from's solid faces by z plane in parallel (1.95 -> 0.62 s per step at 11.6 M cells), the f32 fine export derived from the f64 level the refresh already builds (export_fine 1.04 -> 0.34 s); sub-laps in the mask rebuild; rebuild block 9.4 -> 7.4 s

Co-Authored-By: Claude Fable 5.1 <[email protected]>

crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/impose.rs

@@ -5,6 +5,7 @@
 
 use super::body::Body;
 use super::wall::{FaceKind, Mask};
+use rayon::prelude::*;
 
 /// The imposition band of a cut mask, in cells.
 pub(crate) const IMPOSE_BAND_CELLS: f64 = 4.0;
@@ -52,55 +53,70 @@ impl Mask {
             Some(c) => (&c.d_u, &c.d_v, &c.d_w),
             None => (&[], &[], &[]),
         };
-        for k in 0..nz {
+        // One z plane per task: every solid face is written from the body's
-            for j in 0..ny {
+        // surface velocity alone (no reduction), so the parallel loop is
-                for i in 1..nx {
+        // digit-identical to the serial one (PERF-3 P1-1: 1.95 s of a 9.4 s
-                    let idx = g.uface(k, j, i);
+        // step at 11.6 M cells, serial).
-                    if self.u_kind[idx] == FaceKind::Solid && near(d_u, idx) {
+        u.par_chunks_mut(ny * (nx + 1))
-                        u[idx] = body
+            .take(nz)
-                            .surface_velocity(
+            .enumerate()
-                                i as f64 * dx,
+            .for_each(|(k, plane)| {
-                                (j as f64 + 0.5) * dy,
+                for j in 0..ny {
-                                (k as f64 + 0.5) * dz,
+                    for i in 1..nx {
-                                t,
+                        let idx = g.uface(k, j, i);
-                            )
+                        if self.u_kind[idx] == FaceKind::Solid && near(d_u, idx) {
-                            .0;
+                            plane[j * (nx + 1) + i] = body
+                                .surface_velocity(
+                                    i as f64 * dx,
+                                    (j as f64 + 0.5) * dy,
+                                    (k as f64 + 0.5) * dz,
+                                    t,
+                                )
+                                .0;
+                        }
                     }
                 }
-            }
+            });
-            for j in 1..ny {
+        v.par_chunks_mut((ny + 1) * nx)
-                for i in 0..nx {
+            .take(nz)
-                    let idx = g.vface(k, j, i);
+            .enumerate()
-                    if self.v_kind[idx] == FaceKind::Solid && near(d_v, idx) {
+            .for_each(|(k, plane)| {
-                        v[idx] = body
+                for j in 1..ny {
-                            .surface_velocity(
+                    for i in 0..nx {
-                                (i as f64 + 0.5) * dx,
+                        let idx = g.vface(k, j, i);
-                                j as f64 * dy,
+                        if self.v_kind[idx] == FaceKind::Solid && near(d_v, idx) {
-                                (k as f64 + 0.5) * dz,
+                            plane[j * nx + i] = body
-                                t,
+                                .surface_velocity(
-                            )
+                                    (i as f64 + 0.5) * dx,
-                            .1;
+                                    j as f64 * dy,
+                                    (k as f64 + 0.5) * dz,
+                                    t,
+                                )
+                                .1;
+                        }
                     }
                 }
-            }
+            });
-        }
+        w.par_chunks_mut(ny * nx)
-        for k in 0..=nz {
+            .take(nz + 1)
-            for j in 0..ny {
+            .enumerate()
-                for i in 0..nx {
+            .for_each(|(k, plane)| {
-                    let idx = g.wface(k, j, i);
+                for j in 0..ny {
-                    if self.w_kind[idx] == FaceKind::Solid && near(d_w, idx) {
+                    for i in 0..nx {
-                        w[idx] = body
+                        let idx = g.wface(k, j, i);
-                            .surface_velocity(
+                        if self.w_kind[idx] == FaceKind::Solid && near(d_w, idx) {
-                                (i as f64 + 0.5) * dx,
+                            plane[j * nx + i] = body
-                                (j as f64 + 0.5) * dy,
+                                .surface_velocity(
-                                k as f64 * dz,
+                                    (i as f64 + 0.5) * dx,
-                                t,
+                                    (j as f64 + 0.5) * dy,
-                            )
+                                    k as f64 * dz,
-                            .2;
+                                    t,
+                                )
+                                .2;
+                        }
                     }
                 }
-            }
+            });
-        }
         let u_vals: Vec<f64> = self
             .u_ghosts
             .iter()

crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/poisson/device_cg.rs

@@ -255,7 +255,7 @@ impl DeviceCg {
         let l_key = lap.elapsed();
         // The V-cycle's finest level follows the operator (its coarser
         // levels stay): the fine level alone, no hierarchy build.
-        let fine_level = super::export::export_fine(problem);
+        let fine_level = super::export::export_fine_from(&fine);
         let l_export = lap.elapsed();
         self.vcycle.refresh_fine(&fine_level);
         if profile {

crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/poisson/export.rs

@@ -125,6 +125,40 @@ pub fn export_fine(problem: &Problem) -> LevelExport {
     }
 }
 
+/// As [`export_fine`], from the f64 fine level a caller already built
+/// (PERF-3 P1-1): the f32 level is the same masked coefficients cast, so this
+/// is the identical export without a second clone of the problem and a
+/// second level build.
+pub fn export_fine_from(lv: &Level<f64>) -> LevelExport {
+    let (_, coarse_of) = lv.coarsen();
+    let to_u32 = |v: &[usize]| {
+        v.iter()
+            .map(|&i| if i == usize::MAX { u32::MAX } else { i as u32 })
+            .collect::<Vec<u32>>()
+    };
+    let cast = |v: &[f64]| v.iter().map(|&x| x as f32).collect::<Vec<f32>>();
+    LevelExport {
+        nx: lv.problem.nx,
+        ny: lv.problem.ny,
+        nz: lv.problem.nz,
+        cells: to_u32(&lv.cells),
+        red: to_u32(&lv.red),
+        black: to_u32(&lv.black),
+        top: to_u32(&lv.top),
+        bot: to_u32(&lv.bot),
+        coarse_of: to_u32(&coarse_of),
+        children_ptr: Vec::new(),
+        children_idx: Vec::new(),
+        ae: cast(&lv.ae),
+        aw: cast(&lv.aw),
+        an: cast(&lv.an),
+        as_: cast(&lv.as_),
+        at: cast(&lv.at),
+        ab: cast(&lv.ab),
+        ap: cast(&lv.ap),
+    }
+}
+
 /// `z = M⁻¹ r` by the host f32 V-cycle: the reference a device V-cycle is measured against.
 pub fn vcycle_f32_reference(
     problem: &Problem,

crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/moving.rs

@@ -51,6 +51,7 @@ impl Solver {
             }
         }
         let l_step = lap.elapsed();
+        let sub = std::time::Instant::now();
         new_mask.impose_from(
             body,
             &field.u_old,
@@ -61,6 +62,7 @@ impl Solver {
             &mut field.w,
             t_new,
         );
+        let s_impose = sub.elapsed();
         if new_mask.cut().is_some() {
             let (table, correction) = match &self.mask {
                 Some(old_mask) => new_mask.gcl_flux_table(old_mask, dt),
@@ -69,6 +71,7 @@ impl Solver {
             self.wall_fluxes = table;
             self.last_ghost_correction = correction;
         }
+        let s_gcl = sub.elapsed();
         if let Some(old) = &self.mask {
             self.apertures_old = old
                 .cut()
@@ -86,10 +89,13 @@ impl Solver {
         if std::env::var("RTX_E3_MOVING_PROFILE").is_ok() {
             let ms = |d: std::time::Duration| d.as_secs_f64() * 1e3;
             eprintln!(
-                "    mask laps: build_mask {:.0} ms, refill + step apertures + merging {:.0} ms, impose + GCL + volumes {:.0} ms",
+                "    mask laps: build_mask {:.0} ms, refill + step apertures + merging {:.0} ms, impose + GCL + volumes {:.0} ms (impose {:.0}, GCL table {:.0}, old apertures + volumes {:.0})",
                 ms(l_build),
                 ms(l_step - l_build),
-                ms(lap.elapsed() - l_step)
+                ms(lap.elapsed() - l_step),
+                ms(s_impose),
+                ms(s_gcl - s_impose),
+                ms(sub.elapsed() - s_gcl)
             );
         }
         self.mask = Some(new_mask);