Skip to content

Conversation

@sbryngelson
Copy link
Member

@sbryngelson sbryngelson commented Dec 26, 2025

PR Type

Enhancement, Tests


Description

  • Implement multi-mode perturbation for channel flow initialization with configurable amplitude

  • Add skew-symmetric convective scheme for energy-conserving DNS simulations

  • Enhance configuration system with perturbation amplitude parameter and inline comment parsing

  • Extend test suite with skew-symmetric scheme validation and new turbulent channel case


Diagram Walkthrough

flowchart LR
  A["Config System"] -->|perturbation_amplitude| B["Perturbed Field Init"]
  C["Convective Scheme Enum"] -->|SkewSymmetric| D["Solver Kernels"]
  B -->|multi-mode streamfunction| E["Channel Flow Setup"]
  D -->|energy-conserving form| F["DNS Simulations"]
  G["Adaptive dt"] -->|set_dt method| H["Time Stepping"]
  E -->|unsteady mode| I["Turbulent Case"]
Loading

File Walkthrough

Relevant files
Enhancement
main_channel.cpp
Multi-mode perturbations and improved initialization         

app/main_channel.cpp

  • Refactored create_perturbed_channel_field() to support multi-mode
    perturbations with weighted spectrum (1/k² scaling)
  • Added configurable base flow parameters (u_base, v_base) for flexible
    initialization
  • Enhanced unsteady mode initialization with conditional perturbation
    application based on config.perturbation_amplitude
  • Fixed adaptive dt application by calling solver.set_dt() instead of
    discarding return value
  • Added VTK output diagnostics reporting current dt and CFL number
  • Extended VTK output to include u and v velocity components as separate
    scalar fields
+64/-12 
config.hpp
Configuration enhancements for schemes and perturbations 

include/config.hpp

  • Added SkewSymmetric option to ConvectiveScheme enum with documentation
  • Added perturbation_amplitude configuration parameter (default 0.0) for
    unsteady mode
  • Updated enum documentation for existing convective schemes
+4/-2     
solver.hpp
Add setter for adaptive time step control                               

include/solver.hpp

  • Added set_dt() method to allow external control of current time step
    for adaptive time stepping
+3/-0     
config.cpp
Config parsing for perturbations and skew-symmetric scheme

src/config.cpp

  • Implemented inline comment stripping in config file parser (removes
    text after '#')
  • Added support for parsing perturbation_amplitude parameter from config
    files
  • Extended command-line argument parsing for --perturbation_amplitude /
    --perturb_amp flags
  • Added skew-symmetric scheme recognition in both config file and CLI
    parsing ("skew", "skew_symmetric", "skew-symmetric")
  • Updated help message to document perturbation amplitude option
  • Enhanced config print output to display skew-symmetric scheme name
+21/-1   
solver.cpp
Skew-symmetric convective scheme implementation                   

src/solver.cpp

  • Refactored convective kernels to accept scheme enum (0=Central,
    1=Upwind, 2=SkewSymmetric) instead of boolean
  • Implemented skew-symmetric form for u-momentum: 0.5[(u·∇)u +
    ∂/∂x(½|u|²)] with kinetic energy gradient computation
  • Implemented skew-symmetric form for v-momentum: 0.5[(u·∇)v +
    ∂/∂y(½|u|²)] with kinetic energy gradient computation
  • Updated OpenMP target pragmas to pass scheme parameter to kernels
  • Added CFL and dt diagnostics output when writing VTK files in adaptive
    dt mode
+100/-9 
Tests
test_solver_cpu_gpu.cpp
Comprehensive skew-symmetric scheme test coverage               

tests/test_solver_cpu_gpu.cpp

  • Added Case D in metrics collection for skew-symmetric scheme
    validation on periodic domain
  • Added test_taylor_green_skew_symmetric() test comparing CPU/GPU
    solvers with skew-symmetric scheme
  • Added test_channel_skew_symmetric() test for channel flow with
    skew-symmetric scheme
  • Updated test_various_grids() to use skew-symmetric scheme and renamed
    to Test 5
  • Updated main test runner to execute new skew-symmetric scheme tests
+191/-3 
Configuration changes
turbulent.cfg
New DNS turbulent channel flow configuration                         

examples/channel/cases/turbulent.cfg

  • New DNS-like turbulent channel flow configuration file
  • High-resolution grid (256×512) with wall-normal stretching for
    turbulence resolution
  • Elevated Reynolds number (Re_tau ≈ 300) with strong pressure gradient
  • Skew-symmetric convective scheme for energy conservation
  • Adaptive time stepping with CFL_max=0.3 and small base dt (1e-5)
  • Extended simulation (5M iterations) with moderate perturbation
    amplitude (8.0)
  • Frequent snapshot output (250 snapshots) for turbulence evolution
    tracking
+52/-0   

Summary by CodeRabbit

  • New Features

    • Added Skew-Symmetric convective scheme option and new perturbation_amplitude config/CLI flag
    • Multi-mode perturbation initialization on top of configurable base flow; solver can accept explicit set_dt so adaptive dt is applied
    • VTK output now exports u/v velocities, pressure gradient, and prints dt/CFL diagnostics
    • New DNS-like turbulent channel example config
  • Tests

    • Replaced legacy CPU-GPU parity tests with a two-build dump/compare workflow including skew-symmetric case
  • Chores

    • Added two-build test orchestration script

✏️ Tip: You can customize this high-level summary in your review settings.

@sbryngelson sbryngelson requested a review from Copilot December 26, 2025 21:32
@coderabbitai
Copy link

coderabbitai bot commented Dec 26, 2025

Note

Other AI code review bot(s) detected

CodeRabbit has detected other AI code review bot(s) in this pull request and will avoid duplicating their findings in the review comments. This may lead to a less comprehensive review.

📝 Walkthrough

Walkthrough

Adds SkewSymmetric convective option, multi‑mode divergence‑free perturbation initialization with configurable base flow and amplitude, propagates adaptive dt via a new setter, extends config/CLI parsing for perturbation and scheme aliases, updates solver kernels and VTK output, and adds a two‑build test workflow.

Changes

Cohort / File(s) Summary
Configuration & Public API
include/config.hpp, src/config.cpp, include/solver.hpp
Adds ConvectiveScheme::SkewSymmetric; new perturbation_amplitude config field with CLI alias --perturb_amp; config parser strips inline # comments and validates numeric CLI args; adds RANSSolver::set_dt(double).
Convection Kernel & Solver I/O
src/solver.cpp
Replaces boolean use_central with integer scheme param across convective face kernels and callers; implements skew‑symmetric convective form with per‑face K and derivatives and boundary fallbacks; VTK writer now exports u_velocity/v_velocity, pressure_gradient, and logs dt/CFL when adaptive_dt is active.
Channel Initialization & Time‑stepping
app/main_channel.cpp
Introduces multi‑mode divergence‑free perturbation (N_modes=4, 1/k^2 weighting, deterministic phases); create_perturbed_channel_field now accepts u_base/v_base; unsteady init applies base+perturbation when perturbation_amplitude>0; adaptive dt is applied via solver.set_dt(solver.compute_adaptive_dt()).
Examples / Cases
examples/channel/cases/turbulent.cfg
Adds a DNS‑style turbulent channel config: unsteady, high resolution, SkewSymmetric scheme, adaptive dt/CFL settings, perturbation_amplitude, Poisson tolerances, and output controls.
Tests & CI helpers
tests/test_solver_cpu_gpu.cpp, run_two_build_test.sh
Reworks CPU↔GPU parity tests to a two‑build dump/compare workflow with metric files; adds run_two_build_test.sh to build CPU reference and run GPU comparison; consolidates metric collection and removes legacy per-test comparisons.

Sequence Diagram(s)

sequenceDiagram
  autonumber
  participant User as User/CLI
  participant Config as Config Loader
  participant App as app/main_channel
  participant Mesh as Mesh Builder
  participant Solver as RANSSolver
  participant IO as VTK/Output

  User->>Config: load config (scheme, perturbation_amplitude)
  Config->>App: provide parsed Config
  App->>Mesh: build mesh
  App->>App: create_perturbed_channel_field(mesh, amplitude, u_base, v_base)
  App->>Solver: initialize(base or base+perturbation)
  loop Unsteady steps
    Solver->>Solver: compute_adaptive_dt()
    Solver->>Solver: set_dt(dt)    %% propagate computed dt
    Solver->>Solver: advance timestep (convective scheme selects Central/Upwind/SkewSymmetric)
    Solver->>IO: write_vtk(u_velocity, v_velocity, pressure_gradient, dt/CFL diagnostics)
  end
Loading

Estimated code review effort

🎯 4 (Complex) | ⏱️ ~45 minutes

Poem

🐰 I mulled in four small modes and leapt,

atop the base where soft eddies swept,
skewed currents hummed a careful tune,
timesteps settled by afternoon,
tests hopped home — the patch is kept. 🥕

Pre-merge checks and finishing touches

❌ Failed checks (1 warning, 1 inconclusive)
Check name Status Explanation Resolution
Docstring Coverage ⚠️ Warning Docstring coverage is 45.00% which is insufficient. The required threshold is 80.00%. You can run @coderabbitai generate docstrings to improve docstring coverage.
Title check ❓ Inconclusive The title uses vague, informal language ('multmodal stuff') that lacks clarity and professionalism, making it difficult for reviewers scanning history to understand the specific technical changes. Revise the title to be more descriptive and professional, e.g., 'Add multi-mode perturbations and skew-symmetric convection scheme' or 'Implement configurable multi-mode flow perturbations and energy-conserving DNS support'.
✅ Passed checks (1 passed)
Check name Status Explanation
Description Check ✅ Passed Check skipped - CodeRabbit’s high-level summary is enabled.
✨ Finishing touches
  • 📝 Generate docstrings
🧪 Generate unit tests (beta)
  • Create PR with unit tests
  • Post copyable unit tests in a comment
  • Commit unit tests in branch perturb

Comment @coderabbitai help to get the list of available commands and usage tips.

@qodo-code-review
Copy link

qodo-code-review bot commented Dec 26, 2025

PR Compliance Guide 🔍

Below is a summary of compliance checks for this PR:

Security Compliance
🟢
No security concerns identified No security vulnerabilities detected by AI analysis. Human verification advised for critical code.
Ticket Compliance
🎫 No ticket provided
  • Create ticket/issue
Codebase Duplication Compliance
Codebase context is not defined

Follow the guide to enable codebase context checks.

Custom Compliance
🟢
Generic: Comprehensive Audit Trails

Objective: To create a detailed and reliable record of critical system actions for security analysis
and compliance.

Status: Passed

Learn more about managing compliance generic rules or creating your own custom rules

Generic: Secure Error Handling

Objective: To prevent the leakage of sensitive system information through error messages while
providing sufficient detail for internal debugging.

Status: Passed

Learn more about managing compliance generic rules or creating your own custom rules

🔴
Generic: Meaningful Naming and Self-Documenting Code

Objective: Ensure all identifiers clearly express their purpose and intent, making code
self-documenting

Status:
Magic scheme integers: New code uses raw integer values (scheme == 0/1/2) to represent convective schemes, which
is not self-documenting and obscures intent compared to using the ConvectiveScheme enum.

Referred Code
const bool use_central = (scheme == 0 || scheme == 2);   // 0: Central, 1: Upwind, 2: SkewSymmetric
const bool use_skew    = (scheme == 2);

Learn more about managing compliance generic rules or creating your own custom rules

Generic: Robust Error Handling and Edge Case Management

Objective: Ensure comprehensive error handling that provides meaningful context and graceful
degradation

Status:
Incorrect comment parsing: Inline comment stripping claims to ignore # inside quotes but the implementation always
truncates at the first #, which can misparse valid quoted values and is an unhandled edge
case.

Referred Code
// Strip inline comments (everything after '#' that's not inside quotes)
auto comment_pos = value.find('#');
if (comment_pos != std::string::npos) {
    value = trim(value.substr(0, comment_pos));
}

Learn more about managing compliance generic rules or creating your own custom rules

Generic: Secure Logging Practices

Objective: To ensure logs are useful for debugging and auditing without exposing sensitive
information like PII, PHI, or cardholder data.

Status:
Unstructured console logs: The new VTK diagnostics are written as free-form std::cout strings rather than structured
logs, which may conflict with environments requiring structured (e.g., JSON) logging.

Referred Code
// Diagnostics: report current dt and CFL whenever we write a VTK file
if (config_.adaptive_dt) {
    // u_max based on cell-centered magnitude (host-side; data just synced above in GPU builds)
    double u_max = velocity_.max_magnitude();
    double dx_min = std::min(mesh_->dx, mesh_->dy);
    double cfl = (u_max > 1e-10) ? (u_max * current_dt_ / dx_min) : 0.0;

    std::cout << "[VTK] dt = " << std::scientific << std::setprecision(3) << current_dt_
              << ", CFL = " << std::fixed << std::setprecision(3) << cfl
              << " (CFL_max = " << config_.CFL_max
              << "), file = " << filename << "\n";
}

Learn more about managing compliance generic rules or creating your own custom rules

Generic: Security-First Input Validation and Data Handling

Objective: Ensure all data inputs are validated, sanitized, and handled securely to prevent
vulnerabilities

Status:
CLI value parsing: The new --perturbation_amplitude/--perturb_amp parsing uses std::stod without visible
exception handling or range validation, which may allow crashes or invalid values
depending on surrounding code behavior not shown in the diff.

Referred Code
} else if (arg == "--perturbation_amplitude" || arg == "--perturb_amp") {
    if (i + 1 < argc) {
        perturbation_amplitude = std::stod(argv[++i]);
    } else {
        std::cerr << "ERROR: --perturbation_amplitude requires a value\n";
        return;
    }

Learn more about managing compliance generic rules or creating your own custom rules

  • Update
Compliance status legend 🟢 - Fully Compliant
🟡 - Partial Compliant
🔴 - Not Compliant
⚪ - Requires Further Human Verification
🏷️ - Compliance label

@qodo-code-review
Copy link

qodo-code-review bot commented Dec 26, 2025

PR Code Suggestions ✨

Latest suggestions up to ee38b5c

CategorySuggestion                                                                                                                                    Impact
Possible issue
Correct streamfunction x-derivative
Suggestion Impact:The commit updates the mode summation for dpsi_dx to include a multiplicative kx factor (A_m * kx), matching the intended correction for the x-derivative in the v component, and adds an explanatory comment.

code diff:

@@ -180,7 +180,8 @@
                 double w_m = (1.0 / (m * m)) / weight_sum;
                 double A_m = amplitude * std::sqrt(w_m);
                 double phase = 0.5 * m;  // Deterministic phase offset
-                sum_modes += A_m * std::cos(kx * x + phase);
+                // Include kx factor for ∂/∂x derivative: v = -∂ψ/∂x
+                sum_modes += (A_m * kx) * std::cos(kx * x + phase);

Correct the calculation of the v velocity component from the streamfunction by
multiplying by the wavenumber kx. This fixes a mathematical error that violates
the divergence-free property of the initial velocity field.

app/main_channel.cpp [166-189]

 // v = -dψ/dx at y-faces
 for (int j = mesh.j_begin(); j <= mesh.j_end(); ++j) {
     // Use face coordinate array (correct for stretched grids)
     double y = mesh.yf[j];
     double s = std::sin(0.5 * M_PI * (y + 1.0));
     double s2 = s * s;  // Wall factor: 0 at y=±1
 
     for (int i = mesh.i_begin(); i < mesh.i_end(); ++i) {
         double x = mesh.x(i);
 
         // Sum over multiple modes
         double sum_modes = 0.0;
         for (int m = 1; m <= N_modes; ++m) {
             double kx = m * 2.0 * M_PI / Lx;
             double w_m = (1.0 / (m * m)) / weight_sum;
             double A_m = amplitude * std::sqrt(w_m);
             double phase = 0.5 * m;  // Deterministic phase offset
-            sum_modes += A_m * std::cos(kx * x + phase);
+            sum_modes += (A_m * kx) * std::cos(kx * x + phase);
         }
 
         double dpsi_dx = sum_modes * s2;
         vel.v(i, j) += -dpsi_dx;
     }
 }

[Suggestion processed]

Suggestion importance[1-10]: 9

__

Why: The suggestion correctly identifies a mathematical error in the calculation of the v velocity component from the streamfunction, which breaks the divergence-free property of the initial velocity field. This is a critical bug in the new perturbation logic.

High
  • Update

Previous suggestions

✅ Suggestions up to commit bdd6b0b
CategorySuggestion                                                                                                                                    Impact
Possible issue
Prevent out-of-bounds memory access
Suggestion Impact:The commit addresses the same out-of-bounds risk in the skew-symmetric stencil, but instead of increasing the mesh ghost-layer count to 2, it adds boundary guards (i>=2, j>=2) and falls back to the pure advective form near the west/south boundaries to avoid underflowing indices (i-2, j-2).

code diff:

     if (use_skew) {
-        // Skew-symmetric form: 0.5[(u·∇)u + ∂/∂x(½|u|²)] for u-momentum
-        // East face (i+1)
-        const double uu_e = u_ptr[j * u_stride + (i+1)];
-        const double v_bl_e = v_ptr[j * v_stride + i];
-        const double v_br_e = v_ptr[j * v_stride + (i+1)];
-        const double v_tl_e = v_ptr[(j+1) * v_stride + i];
-        const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
-        const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
-        const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);
-        
-        // West face (i-1)
-        const double uu_w = u_ptr[j * u_stride + (i-1)];
-        const double v_bl_w = v_ptr[j * v_stride + (i-2)];
-        const double v_br_w = v_ptr[j * v_stride + (i-1)];
-        const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
-        const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
-        const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
-        const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
-        
-        const double dKdx = (K_e - K_w) / (2.0 * dx);
-        conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
+        // Skew-symmetric form: 0.5[(u·∇)u + ∂/∂x(½|u|²)] for u-momentum.
+        // This stencil accesses v at (i-2), so with Nghost = 1 the first
+        // interior face (i == 1) would underflow. To avoid out-of-bounds
+        // access we fall back to the advective form for i < 2.
+        if (i >= 2) {
+            // East face (i+1)
+            const double uu_e = u_ptr[j * u_stride + (i+1)];
+            const double v_bl_e = v_ptr[j * v_stride + i];
+            const double v_br_e = v_ptr[j * v_stride + (i+1)];
+            const double v_tl_e = v_ptr[(j+1) * v_stride + i];
+            const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
+            const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
+            const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);
+            
+            // West face (i-1)
+            const double uu_w = u_ptr[j * u_stride + (i-1)];
+            const double v_bl_w = v_ptr[j * v_stride + (i-2)];
+            const double v_br_w = v_ptr[j * v_stride + (i-1)];
+            const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
+            const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
+            const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
+            const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
+            
+            const double dKdx = (K_e - K_w) / (2.0 * dx);
+            conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
+        } else {
+            // Near the west boundary: use pure advective form to stay in-bounds.
+            conv_u_ptr[conv_idx] = adv_u;
+        }
     } else {
         conv_u_ptr[conv_idx] = adv_u;
     }
@@ -493,29 +501,37 @@
     const double adv_v = uu * dvdx + vv * dvdy;
     
     if (use_skew) {
-        // Skew-symmetric form: 0.5[(u·∇)v + ∂/∂y(½|u|²)] for v-momentum
-        // North face (j+1)
-        const int v_idx_n = (j+1) * v_stride + i;
-        const double vv_n = v_ptr[v_idx_n];
-        const double u_bl_n = u_ptr[j * u_stride + i];
-        const double u_br_n = u_ptr[j * u_stride + (i+1)];
-        const double u_tl_n = u_ptr[(j+1) * u_stride + i];
-        const double u_tr_n = u_ptr[(j+1) * u_stride + (i+1)];
-        const double uu_n = 0.25 * (u_bl_n + u_br_n + u_tl_n + u_tr_n);
-        const double K_n = 0.5 * (uu_n * uu_n + vv_n * vv_n);
-        
-        // South face (j-1)
-        const int v_idx_s = (j-1) * v_stride + i;
-        const double vv_s = v_ptr[v_idx_s];
-        const double u_bl_s = u_ptr[(j-2) * u_stride + i];
-        const double u_br_s = u_ptr[(j-2) * u_stride + (i+1)];
-        const double u_tl_s = u_ptr[(j-1) * u_stride + i];
-        const double u_tr_s = u_ptr[(j-1) * u_stride + (i+1)];
-        const double uu_s = 0.25 * (u_bl_s + u_br_s + u_tl_s + u_tr_s);
-        const double K_s = 0.5 * (uu_s * uu_s + vv_s * vv_s);
-        
-        const double dKdy = (K_n - K_s) / (2.0 * dy);
-        conv_v_ptr[conv_idx] = 0.5 * (adv_v + dKdy);
+        // Skew-symmetric form: 0.5[(u·∇)v + ∂/∂y(½|u|²)] for v-momentum.
+        // This stencil accesses u at (j-2), so with Nghost = 1 the first
+        // interior row (j == 1) would underflow. To avoid out-of-bounds
+        // access we fall back to the advective form for j < 2.
+        if (j >= 2) {
+            // North face (j+1)
+            const int v_idx_n = (j+1) * v_stride + i;
+            const double vv_n = v_ptr[v_idx_n];
+            const double u_bl_n = u_ptr[j * u_stride + i];
+            const double u_br_n = u_ptr[j * u_stride + (i+1)];
+            const double u_tl_n = u_ptr[(j+1) * u_stride + i];
+            const double u_tr_n = u_ptr[(j+1) * u_stride + (i+1)];
+            const double uu_n = 0.25 * (u_bl_n + u_br_n + u_tl_n + u_tr_n);
+            const double K_n = 0.5 * (uu_n * uu_n + vv_n * vv_n);
+            
+            // South face (j-1)
+            const int v_idx_s = (j-1) * v_stride + i;
+            const double vv_s = v_ptr[v_idx_s];
+            const double u_bl_s = u_ptr[(j-2) * u_stride + i];
+            const double u_br_s = u_ptr[(j-2) * u_stride + (i+1)];
+            const double u_tl_s = u_ptr[(j-1) * u_stride + i];
+            const double u_tr_s = u_ptr[(j-1) * u_stride + (i+1)];
+            const double uu_s = 0.25 * (u_bl_s + u_br_s + u_tl_s + u_tr_s);
+            const double K_s = 0.5 * (uu_s * uu_s + vv_s * vv_s);
+            
+            const double dKdy = (K_n - K_s) / (2.0 * dy);
+            conv_v_ptr[conv_idx] = 0.5 * (adv_v + dKdy);
+        } else {
+            // Near the south boundary: use pure advective form to stay in-bounds.
+            conv_v_ptr[conv_idx] = adv_v;
+        }

Prevent an out-of-bounds memory access in the skew-symmetric convection kernel
by increasing the number of ghost cell layers to at least 2 when initializing
the mesh for this scheme.

src/solver.cpp [433-440]

-// West face (i-1)
-const double uu_w = u_ptr[j * u_stride + (i-1)];
-const double v_bl_w = v_ptr[j * v_stride + (i-2)];
-const double v_br_w = v_ptr[j * v_stride + (i-1)];
-const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
-const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
-const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
-const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
+// This suggestion is conceptual and requires a change in a different location.
+// The fix is to ensure the mesh is created with at least 2 ghost layers
+// when using the skew-symmetric scheme. For example, in `mesh.init_uniform()`:
+//
+// int n_ghost = (config.convective_scheme == ConvectiveScheme::SkewSymmetric) ? 2 : 1;
+// mesh.init_uniform(..., n_ghost);
+//
+// The code snippet below remains unchanged, but the underlying mesh data it accesses
+// must have a larger ghost cell region to prevent out-of-bounds reads.
 
+    // West face (i-1)
+    const double uu_w = u_ptr[j * u_stride + (i-1)];
+    const double v_bl_w = v_ptr[j * v_stride + (i-2)];
+    const double v_br_w = v_ptr[j * v_stride + (i-1)];
+    const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
+    const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
+    const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
+    const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
+

[Suggestion processed]

Suggestion importance[1-10]: 10

__

Why: This suggestion correctly identifies a critical out-of-bounds memory access bug in the new skew-symmetric scheme implementation due to an insufficient number of ghost cells. This can lead to crashes or incorrect results and is a major correctness issue.

High
Fix incorrect v-velocity perturbation calculation
Suggestion Impact:The commit updates the summation to multiply A_m by kx (with an explanatory comment), matching the intended chain-rule correction for computing dpsi_dx used in the v-velocity component.

code diff:

-                sum_modes += A_m * std::cos(kx * x + phase);
+                // Include kx factor for ∂/∂x derivative: v = -∂ψ/∂x
+                sum_modes += (A_m * kx) * std::cos(kx * x + phase);

Fix the calculation of the v velocity component in
create_perturbed_channel_field by including the missing kx factor in the
summation, ensuring the resulting velocity field is divergence-free.

app/main_channel.cpp [176-187]

 // Sum over multiple modes
 double sum_modes = 0.0;
 for (int m = 1; m <= N_modes; ++m) {
     double kx = m * 2.0 * M_PI / Lx;
     double w_m = (1.0 / (m * m)) / weight_sum;
     double A_m = amplitude * std::sqrt(w_m);
     double phase = 0.5 * m;  // Deterministic phase offset
-    sum_modes += A_m * std::cos(kx * x + phase);
+    sum_modes += A_m * kx * std::cos(kx * x + phase);
 }
 
 double dpsi_dx = sum_modes * s2;
 vel.v(i, j) += -dpsi_dx;
Suggestion importance[1-10]: 9

__

Why: The suggestion correctly identifies a mathematical error in the implementation of v = -∂ψ/∂x, where the kx factor from the chain rule is missing. This bug breaks the divergence-free property of the initial velocity field, which is critical for the solver's correctness.

High
High-level
Refactor the skew-symmetric scheme implementation

The implemented skew-symmetric scheme is a non-standard form that may not
conserve energy. It should be replaced with the standard formulation, which
averages the advective and divergence forms (0.5 * [(u⋅∇)φ + ∇⋅(uφ)]), to ensure
discrete energy conservation for DNS.

Examples:

src/solver.cpp [422-446]
    if (use_skew) {
        // Skew-symmetric form: 0.5[(u·∇)u + ∂/∂x(½|u|²)] for u-momentum
        // East face (i+1)
        const double uu_e = u_ptr[j * u_stride + (i+1)];
        const double v_bl_e = v_ptr[j * v_stride + i];
        const double v_br_e = v_ptr[j * v_stride + (i+1)];
        const double v_tl_e = v_ptr[(j+1) * v_stride + i];
        const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
        const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
        const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);

 ... (clipped 15 lines)
src/solver.cpp [495-521]
    if (use_skew) {
        // Skew-symmetric form: 0.5[(u·∇)v + ∂/∂y(½|u|²)] for v-momentum
        // North face (j+1)
        const int v_idx_n = (j+1) * v_stride + i;
        const double vv_n = v_ptr[v_idx_n];
        const double u_bl_n = u_ptr[j * u_stride + i];
        const double u_br_n = u_ptr[j * u_stride + (i+1)];
        const double u_tl_n = u_ptr[(j+1) * u_stride + i];
        const double u_tr_n = u_ptr[(j+1) * u_stride + (i+1)];
        const double uu_n = 0.25 * (u_bl_n + u_br_n + u_tl_n + u_tr_n);

 ... (clipped 17 lines)

Solution Walkthrough:

Before:

// In convective_u_face_kernel_staggered
if (use_skew) {
    // Implements: 0.5 * [(u·∇)u + ∂/∂x(½|u|²)]
    const double adv_u = ...; // (u·∇)u
    const double K_e = 0.5 * (uu_e*uu_e + vv_e*vv_e);
    const double K_w = 0.5 * (uu_w*uu_w + vv_w*vv_w);
    const double dKdx = (K_e - K_w) / (2.0 * dx);
    conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
}

// In convective_v_face_kernel_staggered
if (use_skew) {
    // Implements: 0.5 * [(u·∇)v + ∂/∂y(½|u|²)]
    const double adv_v = ...; // (u·∇)v
    const double dKdy = ...;
    conv_v_ptr[conv_idx] = 0.5 * (adv_v + dKdy);
}

After:

// In convective_u_face_kernel_staggered
if (use_skew) {
    // Implements: 0.5 * [(u·∇)u + ∇·(uu)]
    const double adv_u = ...; // (u·∇)u
    const double flux_u_e = u_e * u_e;
    const double flux_u_w = u_w * u_w;
    const double flux_v_n = v_n * u_n;
    const double flux_v_s = v_s * u_s;
    const double div_form = (flux_u_e - flux_u_w)/dx + (flux_v_n - flux_v_s)/dy;
    conv_u_ptr[conv_idx] = 0.5 * (adv_u + div_form);
}

// In convective_v_face_kernel_staggered
if (use_skew) {
    // Implements: 0.5 * [(u·∇)v + ∇·(uv)]
    const double adv_v = ...; // (u·∇)v
    const double div_form = ...; // ∇·(uv)
    conv_v_ptr[conv_idx] = 0.5 * (adv_v + div_form);
}
Suggestion importance[1-10]: 9

__

Why: The suggestion correctly identifies that the implemented skew-symmetric scheme is a non-standard form and may not be discretely energy-conserving, which undermines a primary goal of the PR.

High

Copy link

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 3

🧹 Nitpick comments (3)
src/config.cpp (2)

56-60: Inline comment stripping doesn't handle quoted strings as documented.

The comment mentions "everything after '#' that's not inside quotes" but the implementation doesn't actually check for quotes. This could incorrectly strip valid content from string values containing #.

🔎 Example problematic case

If a config line is:

nn_weights_path = /path/to/#model  # comment

The value would be incorrectly truncated to /path/to/.

Consider whether this edge case matters for your config values. If string values never contain #, the current implementation is fine.


280-286: Consider using std::exit(1) instead of return for consistency with other error paths.

Other error handling in this function (e.g., the --help path at line 324) uses std::exit(0). The return here silently continues parsing remaining arguments, which may lead to unexpected behavior if the user misses the error message.

🔎 Suggested fix
         } else if (arg == "--perturbation_amplitude" || arg == "--perturb_amp") {
             if (i + 1 < argc) {
                 perturbation_amplitude = std::stod(argv[++i]);
             } else {
                 std::cerr << "ERROR: --perturbation_amplitude requires a value\n";
-                return;
+                std::exit(1);
             }
tests/test_solver_cpu_gpu.cpp (1)

673-700: Note: Grid sweep test now always uses skew-symmetric scheme.

This test previously used the default Central scheme and now uses SkewSymmetric. This changes the test's purpose from general grid-size validation to skew-symmetric-specific validation.

Consider whether you want to retain a Central scheme grid sweep test for regression coverage, or if the existing Taylor-Green and channel tests (which use Central) provide sufficient coverage.

📜 Review details

Configuration used: defaults

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between adbe039 and bdd6b0b.

📒 Files selected for processing (7)
  • app/main_channel.cpp
  • examples/channel/cases/turbulent.cfg
  • include/config.hpp
  • include/solver.hpp
  • src/config.cpp
  • src/solver.cpp
  • tests/test_solver_cpu_gpu.cpp
🧰 Additional context used
📓 Path-based instructions (3)
**/*.{cpp,cc,cxx,h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx,h,hpp}: Avoid overly strict floating-point comparisons (e.g., == for doubles); use appropriate tolerances based on algorithm
Fix all compiler warnings before pushing; use [[maybe_unused]] for intentionally unused variables in assertions
Use const for read-only references and mark methods const if they don't modify state
Use RAII for resource management; use smart pointers or RAII wrappers instead of manual new/delete or malloc/free
Use OmpDeviceBuffer wrapper for GPU buffer management instead of manual allocation
Check return values and use exceptions for error conditions with informative error messages

Files:

  • include/solver.hpp
  • app/main_channel.cpp
  • tests/test_solver_cpu_gpu.cpp
  • include/config.hpp
  • src/solver.cpp
  • src/config.cpp
include/**/*.{h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

Every public function needs documentation explaining parameters, return values, and side effects; include usage examples for complex functions

Files:

  • include/solver.hpp
  • include/config.hpp
**/*.{cpp,cc,cxx}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx}: Do not add platform-specific tolerances with #ifdef __APPLE__ or similar platform checks in tests
Include code comments explaining WHY, not WHAT; document numerical algorithms; note non-obvious optimizations

Files:

  • app/main_channel.cpp
  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • src/config.cpp
🧠 Learnings (6)
📓 Common learnings
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability

Applied to files:

  • app/main_channel.cpp
  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • examples/channel/cases/turbulent.cfg
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Always call `set_body_force()` for driven flows; initialize velocity field before solving; set turbulence model before first step

Applied to files:

  • app/main_channel.cpp
  • examples/channel/cases/turbulent.cfg
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Use `solve_steady_with_snapshots()` for automatic VTK output; specify `num_snapshots` in config; files are numbered sequentially plus final

Applied to files:

  • app/main_channel.cpp
  • src/solver.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_tg_validation.cpp : Taylor-Green Vortex test must validate energy decay with <5% error over 100 timesteps using viscous terms, time integration, and periodic BCs

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: For GPU-related changes, also run `./test_before_ci_gpu.sh` to run the complete GPU CI test suite locally

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
🧬 Code graph analysis (3)
app/main_channel.cpp (2)
include/mesh.hpp (4)
  • j (74-74)
  • j (74-74)
  • j (78-78)
  • j (78-78)
include/fields.hpp (12)
  • i (19-19)
  • i (19-19)
  • i (20-20)
  • i (20-20)
  • i (62-62)
  • i (62-62)
  • i (63-63)
  • i (63-63)
  • i (67-67)
  • i (67-67)
  • i (68-68)
  • i (68-68)
src/solver.cpp (3)
include/fields.hpp (18)
  • i (19-19)
  • i (19-19)
  • i (20-20)
  • i (20-20)
  • i (62-62)
  • i (62-62)
  • i (63-63)
  • i (63-63)
  • i (67-67)
  • i (67-67)
  • i (68-68)
  • i (68-68)
  • i (86-86)
  • i (95-95)
  • i (96-96)
  • i (116-116)
  • filename (46-46)
  • filename (101-101)
include/solver.hpp (3)
  • velocity_ (105-105)
  • velocity_ (112-112)
  • filename (137-137)
app/main_periodic_hills.cpp (3)
  • filename (277-279)
  • filename (277-277)
  • file (294-294)
src/config.cpp (1)
include/config.hpp (1)
  • argc (113-113)
⏰ Context from checks skipped due to timeout of 90000ms. You can increase the timeout in your CodeRabbit configuration to a maximum of 15 minutes (900000ms). (3)
  • GitHub Check: Agent
  • GitHub Check: build-and-test (ubuntu-latest, Debug)
  • GitHub Check: gpu-tests
🔇 Additional comments (16)
include/config.hpp (2)

27-32: LGTM: Well-documented convective scheme enum extension.

The enum is clearly documented with each scheme's characteristics (advective form, robustness, energy-conserving). The naming is consistent with CFD conventions.


69-69: LGTM: Clear perturbation amplitude parameter.

The new configuration field is well-documented and has a sensible default (0.0 = no perturbation).

include/solver.hpp (1)

101-102: LGTM: Clean setter for adaptive time stepping.

The setter enables explicit control of the time step, which is necessary for the adaptive dt workflow in the main loop. The implementation is straightforward.

src/config.cpp (1)

122-122: LGTM: Consistent handling of new configuration options.

The perturbation_amplitude and SkewSymmetric scheme are consistently handled across file loading, argument parsing, and status display.

Also applies to: 128-129, 268-269, 472-473

examples/channel/cases/turbulent.cfg (2)

52-52: Confirm perturbation amplitude of 8.0 is appropriate.

An amplitude of 8.0 is quite large compared to typical initialization perturbations (often O(0.01-0.1)). For DNS at Re_tau~300 where velocity scales are O(1), this could introduce unrealistically large initial disturbances.

However, if this is intentional to quickly trigger transition, it may be acceptable. Please verify the physical rationale.


1-52: Useful DNS configuration example.

The configuration demonstrates proper usage of the new features: unsteady mode, skew-symmetric scheme, and perturbation initialization. The comments are helpful for understanding the parameter choices.

src/solver.cpp (4)

495-521: Same out-of-bounds concern for skew-symmetric v-momentum kernel.

The south-face velocity interpolation at lines 510-513 accesses u_ptr[(j-2) * u_stride + i] and similar indices. For faces at j = Ng, this accesses j-2 = Ng-2, which could be out of bounds if Ng < 2.


394-396: LGTM: Clean scheme selection logic.

The use of scheme == 0 || scheme == 2 for central differencing (both Central and SkewSymmetric use central stencils) and scheme == 2 for the skew-symmetric adjustment is clear and correct.

Also applies to: 460-462


2457-2468: LGTM: Helpful diagnostic output for adaptive time stepping.

The VTK write now reports current dt and CFL when adaptive_dt is enabled, which aids debugging and monitoring of unsteady simulations.


2500-2517: LGTM: Additional scalar fields in VTK output.

Exporting u_velocity and v_velocity as separate scalar fields provides more flexibility for post-processing and visualization.

tests/test_solver_cpu_gpu.cpp (3)

288-346: LGTM: Good addition of skew-symmetric test case to metrics collection.

The Case D test mirrors the Taylor-Green setup from Case A but with ConvectiveScheme::SkewSymmetric, providing a reference for CPU/GPU comparison of the new scheme.


548-616: LGTM: Well-structured skew-symmetric Taylor-Green test.

The test properly initializes both CPU and GPU solvers with identical configurations and compares results after 10 time steps. This validates CPU/GPU consistency for the new convective scheme.


833-835: LGTM: New tests integrated into test flow.

The skew-symmetric tests are properly called in the GPU-enabled test path.

app/main_channel.cpp (3)

118-192: LGTM: Well-designed multi-mode perturbation generator.

The refactored create_perturbed_channel_field function:

  1. Uses a 1/k² energy spectrum (lines 136-139), which is physically motivated for developed turbulence
  2. Maintains divergence-free constraint via streamfunction approach
  3. Supports base flow initialization via new u_base/v_base parameters
  4. Uses deterministic phase offsets for reproducibility

The v-component correctly computes -∂ψ/∂x (line 187), ensuring the divergence-free constraint.


310-326: LGTM: Clean conditional initialization for unsteady mode.

The initialization path correctly branches based on perturbation_amplitude:

  • With perturbation: creates multi-mode field with base flow
  • Without: initializes uniform flow

This follows the learnings: "initialize velocity field before solving".


340-342: Good fix: Now properly propagates adaptive dt to solver.

The change from (void)solver.compute_adaptive_dt() to solver.set_dt(solver.compute_adaptive_dt()) ensures the computed adaptive time step is actually used by the solver. This was likely a bug in the previous implementation.

Comment on lines 47 to 49
# Poisson solver
poisson_tol = 1e-4
poisson_max_iter = 3
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🔴 Critical

🧩 Analysis chain

🏁 Script executed:

#!/bin/bash
# Find all poisson_max_iter settings in the codebase
echo "=== poisson_max_iter in .cfg files ==="
find . -type f -name "*.cfg" -exec grep -n "poisson_max_iter" {} + 2>/dev/null

echo -e "\n=== poisson_max_iter in all files (context) ==="
rg -n "poisson_max_iter" -A 1 -B 1 2>/dev/null | head -50

Repository: comp-physics/cfd-nn

Length of output: 4576


🏁 Script executed:

#!/bin/bash
# Check if there's documentation or code about the Poisson solver and convergence
echo "=== Poisson solver related comments/docs ==="
rg -n "Poisson|poisson" --type=cpp --type=h -C 1 examples/ | grep -i "iter\|converge\|cycle\|multigrid" | head -30

Repository: comp-physics/cfd-nn

Length of output: 107


🏁 Script executed:

#!/bin/bash
# Check the solver implementation to understand what max_iter represents
fd -n "poisson|solver" -t f --type cpp . | head -20
rg -n "max_iter.*poisson|poisson.*max_iter" --type=cpp 2>/dev/null

Repository: comp-physics/cfd-nn

Length of output: 1677


Correct poisson_max_iter = 3 to align with other production configurations.

This value is anomalously low. All other production configs use 5,000–20,000 iterations:

  • Most channel/RANS configs: 10,000
  • Fine grid cases: 20,000
  • Simpler cases: 5,000

The only use of 3 iterations in the codebase is in app/profile_poisson_warmstart.cpp, a performance profiling utility. For a 256×512 grid DNS turbulent channel, inadequate Poisson solver convergence will violate the divergence-free constraint and compromise numerical stability. Use at least 10,000 iterations or justify this as an intentional performance study with clear documentation.

🤖 Prompt for AI Agents
In examples/channel/cases/turbulent.cfg around lines 47 to 49, poisson_max_iter
is set to 3 which is far too low for production turbulent channel runs and will
break divergence-free convergence; change poisson_max_iter to at least 10000 to
match other channel/RANS configs (or to 20000 for fine-grid cases) unless this
file is explicitly intended as a profiling/performance study—if it is a
deliberate performance test, add a comment above the setting documenting that
intent and link to app/profile_poisson_warmstart.cpp so reviewers understand it
is not a production configuration.

Copy link

Copilot AI left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Pull request overview

This PR adds support for a skew-symmetric convective scheme and multi-mode perturbations for initial flow fields in unsteady simulations. The skew-symmetric scheme is energy-conserving and better suited for DNS-like turbulent channel flows compared to standard central differences.

Key changes:

  • Implements skew-symmetric convective scheme (energy-conserving form) as an alternative to central/upwind discretization
  • Adds multi-mode streamfunction-based perturbation generator for triggering turbulent transition
  • Adds configuration parameter perturbation_amplitude to control initial perturbation strength
  • Adds diagnostic output (dt, CFL) to VTK writes and separate u/v velocity scalar fields

Reviewed changes

Copilot reviewed 7 out of 7 changed files in this pull request and generated 8 comments.

Show a summary per file
File Description
include/config.hpp Adds SkewSymmetric enum value and perturbation_amplitude config parameter
src/config.cpp Parses skew-symmetric scheme option and perturbation amplitude; adds inline comment stripping
include/solver.hpp Adds set_dt() method for adaptive time stepping
src/solver.cpp Implements skew-symmetric convective kernel with kinetic energy gradient term; adds VTK diagnostics and velocity component outputs
app/main_channel.cpp Replaces single-mode with multi-mode perturbation using 1/k² spectrum; fixes adaptive dt usage
tests/test_solver_cpu_gpu.cpp Adds three new tests (Case D, Test 3, Test 4) for skew-symmetric scheme validation; updates Test 5 numbering
examples/channel/cases/turbulent.cfg New DNS-like turbulent channel configuration using skew-symmetric scheme with perturbations

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

src/config.cpp Outdated
<< " --stretch Use stretched mesh in y\n"
<< " --adaptive_dt Enable adaptive time stepping\n"
<< " --CFL VALUE Max CFL number for adaptive dt (default 0.5)\n"
<< " --scheme SCHEME Convective scheme: central (default), upwind\n"
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The help text for the --scheme option is incomplete and doesn't mention the new skew_symmetric scheme option. It should include all three valid options: central (default), upwind, and skew_symmetric (or skew-symmetric).

Suggested change
<< " --scheme SCHEME Convective scheme: central (default), upwind\n"
<< " --scheme SCHEME Convective scheme: central (default), upwind, skew_symmetric\n"

Copilot uses AI. Check for mistakes.
Comment on lines +154 to +160
for (int m = 1; m <= N_modes; ++m) {
double kx = m * 2.0 * M_PI / Lx;
double w_m = (1.0 / (m * m)) / weight_sum;
double A_m = amplitude * std::sqrt(w_m);
double phase = 0.5 * m; // Deterministic phase offset
sum_modes += A_m * std::sin(kx * x + phase);
}
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The multi-mode computation repeatedly calculates the same values (kx, w_m, A_m, phase) inside the inner loop across all grid points. Consider hoisting these calculations outside the loops to improve performance, especially for large grids.

Copilot uses AI. Check for mistakes.
src/config.cpp Outdated
Comment on lines 56 to 60
// Strip inline comments (everything after '#' that's not inside quotes)
auto comment_pos = value.find('#');
if (comment_pos != std::string::npos) {
value = trim(value.substr(0, comment_pos));
}
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The inline comment stripping logic doesn't properly handle '#' characters inside quoted strings. For example, if a config value contains '#' within quotes (like a filename with '#'), it would incorrectly truncate the value. Consider implementing proper quote-aware parsing to handle this edge case.

Copilot uses AI. Check for mistakes.
src/solver.cpp Outdated
Comment on lines 424 to 443
// East face (i+1)
const double uu_e = u_ptr[j * u_stride + (i+1)];
const double v_bl_e = v_ptr[j * v_stride + i];
const double v_br_e = v_ptr[j * v_stride + (i+1)];
const double v_tl_e = v_ptr[(j+1) * v_stride + i];
const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);

// West face (i-1)
const double uu_w = u_ptr[j * u_stride + (i-1)];
const double v_bl_w = v_ptr[j * v_stride + (i-2)];
const double v_br_w = v_ptr[j * v_stride + (i-1)];
const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);

const double dKdx = (K_e - K_w) / (2.0 * dx);
conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The skew-symmetric scheme accesses v_ptr at index (i-2) when computing the west face kinetic energy. With the default Nghost=1, this would access index -1 when processing the first interior face (i=Ng=1), which is out of bounds. The skew-symmetric scheme requires at least 2 ghost cell layers (Nghost >= 2) to work correctly, but this requirement is not enforced or documented.

Suggested change
// East face (i+1)
const double uu_e = u_ptr[j * u_stride + (i+1)];
const double v_bl_e = v_ptr[j * v_stride + i];
const double v_br_e = v_ptr[j * v_stride + (i+1)];
const double v_tl_e = v_ptr[(j+1) * v_stride + i];
const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);
// West face (i-1)
const double uu_w = u_ptr[j * u_stride + (i-1)];
const double v_bl_w = v_ptr[j * v_stride + (i-2)];
const double v_br_w = v_ptr[j * v_stride + (i-1)];
const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
const double dKdx = (K_e - K_w) / (2.0 * dx);
conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
//
// NOTE: The west-face stencil accesses v_ptr at (i-2). To avoid
// out-of-bounds accesses when only a single ghost cell is present
// (e.g., i == 1), we only apply the skew-symmetric correction when
// i is sufficiently far from the west boundary. Otherwise, we fall
// back to the advective term computed above.
if (i >= 2) {
// East face (i+1)
const double uu_e = u_ptr[j * u_stride + (i+1)];
const double v_bl_e = v_ptr[j * v_stride + i];
const double v_br_e = v_ptr[j * v_stride + (i+1)];
const double v_tl_e = v_ptr[(j+1) * v_stride + i];
const double v_tr_e = v_ptr[(j+1) * v_stride + (i+1)];
const double vv_e = 0.25 * (v_bl_e + v_br_e + v_tl_e + v_tr_e);
const double K_e = 0.5 * (uu_e * uu_e + vv_e * vv_e);
// West face (i-1)
const double uu_w = u_ptr[j * u_stride + (i-1)];
const double v_bl_w = v_ptr[j * v_stride + (i-2)];
const double v_br_w = v_ptr[j * v_stride + (i-1)];
const double v_tl_w = v_ptr[(j+1) * v_stride + (i-2)];
const double v_tr_w = v_ptr[(j+1) * v_stride + (i-1)];
const double vv_w = 0.25 * (v_bl_w + v_br_w + v_tl_w + v_tr_w);
const double K_w = 0.5 * (uu_w * uu_w + vv_w * vv_w);
const double dKdx = (K_e - K_w) / (2.0 * dx);
conv_u_ptr[conv_idx] = 0.5 * (adv_u + dKdx);
} else {
// Near the west boundary: use the non-skew advective term to
// avoid accessing v_ptr with a negative index.
conv_u_ptr[conv_idx] = adv_u;
}

Copilot uses AI. Check for mistakes.
src/solver.cpp Outdated
Comment on lines 510 to 515
const double u_bl_s = u_ptr[(j-2) * u_stride + i];
const double u_br_s = u_ptr[(j-2) * u_stride + (i+1)];
const double u_tl_s = u_ptr[(j-1) * u_stride + i];
const double u_tr_s = u_ptr[(j-1) * u_stride + (i+1)];
const double uu_s = 0.25 * (u_bl_s + u_br_s + u_tl_s + u_tr_s);
const double K_s = 0.5 * (uu_s * uu_s + vv_s * vv_s);
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The skew-symmetric scheme accesses u_ptr at index (j-2) when computing the south face kinetic energy. With the default Nghost=1, this would access index -1 when processing the first interior face (j=Ng=1), which is out of bounds. The skew-symmetric scheme requires at least 2 ghost cell layers (Nghost >= 2) to work correctly, but this requirement is not enforced or documented.

Copilot uses AI. Check for mistakes.
Comment on lines +394 to +395
const bool use_central = (scheme == 0 || scheme == 2); // 0: Central, 1: Upwind, 2: SkewSymmetric
const bool use_skew = (scheme == 2);
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The scheme enum values (0, 1, 2) are hardcoded in comments and comparisons. Consider using named constants or the enum values directly to improve maintainability and avoid potential errors if the enum order changes.

Copilot uses AI. Check for mistakes.
Comment on lines +460 to +461
const bool use_central = (scheme == 0 || scheme == 2); // 0: Central, 1: Upwind, 2: SkewSymmetric
const bool use_skew = (scheme == 2);
Copy link

Copilot AI Dec 26, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The scheme enum values (0, 1, 2) are hardcoded in comments and comparisons. Consider using named constants or the enum values directly to improve maintainability and avoid potential errors if the enum order changes.

Copilot uses AI. Check for mistakes.
@sbryngelson
Copy link
Member Author

/improve

@qodo-code-review
Copy link

Persistent suggestions updated to latest commit ee38b5c

Copy link

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 1

♻️ Duplicate comments (2)
src/config.cpp (2)

56-60: Inline comment stripping doesn't handle '#' within quoted strings.

The current implementation strips everything after the first '#' character, which would incorrectly truncate values containing '#' inside quotes (e.g., path = "file#name.dat" # comment becomes path = "file).

For typical config files this is rarely an issue, but consider quote-aware parsing if config values may legitimately contain '#'.


315-315: Help text for --scheme doesn't mention skew_symmetric option.

The help text shows --scheme SCHEME Convective scheme: central (default), upwind but doesn't document the new skew_symmetric option.

🔎 Suggested fix
-                      << "  --scheme SCHEME   Convective scheme: central (default), upwind\n"
+                      << "  --scheme SCHEME   Convective scheme: central (default), upwind, skew_symmetric\n"
🧹 Nitpick comments (3)
src/config.cpp (2)

128-129: Consider extracting scheme string parsing into a helper function.

The logic to map scheme strings ("skew", "skew_symmetric", "skew-symmetric") to the enum is duplicated between config file loading (lines 128-129) and CLI parsing (lines 268-269). Extracting this into a helper function would reduce duplication and ensure consistency.

🔎 Proposed refactor

Add a helper function at the top of config.cpp:

static ConvectiveScheme parse_convective_scheme(const std::string& scheme_str) {
    if (scheme_str == "upwind") {
        return ConvectiveScheme::Upwind;
    } else if (scheme_str == "skew" || scheme_str == "skew_symmetric" || scheme_str == "skew-symmetric") {
        return ConvectiveScheme::SkewSymmetric;
    } else {
        return ConvectiveScheme::Central;
    }
}

Then use it in both locations:

     auto scheme_str = get_string("convective_scheme", "central");
-    if (scheme_str == "upwind") {
-        convective_scheme = ConvectiveScheme::Upwind;
-    } else if (scheme_str == "skew" || scheme_str == "skew_symmetric" || scheme_str == "skew-symmetric") {
-        convective_scheme = ConvectiveScheme::SkewSymmetric;
-    } else {
-        convective_scheme = ConvectiveScheme::Central;
-    }
+    convective_scheme = parse_convective_scheme(scheme_str);

Also applies to: 268-269


472-473: Consider simplifying nested ternary operator for readability.

The nested ternary operator for printing the convective scheme is correct but can be harder to read. Consider using a switch statement or helper function for clarity.

🔎 Proposed refactor
-    std::cout << "Convective scheme: " 
-              << (convective_scheme == ConvectiveScheme::Upwind ? "Upwind" : 
-                  (convective_scheme == ConvectiveScheme::SkewSymmetric ? "Skew-Symmetric" : "Central")) << "\n"
+    std::cout << "Convective scheme: ";
+    switch (convective_scheme) {
+        case ConvectiveScheme::Upwind: std::cout << "Upwind"; break;
+        case ConvectiveScheme::SkewSymmetric: std::cout << "Skew-Symmetric"; break;
+        default: std::cout << "Central"; break;
+    }
+    std::cout << "\n"
src/solver.cpp (1)

394-395: Hardcoded scheme enum values reduce maintainability.

The code uses magic numbers (0, 1, 2) to represent convective schemes with inline comments explaining the mapping. This is fragile if the enum definition in config.hpp changes order. Consider passing the enum directly or using named constants.

🔎 Proposed fix (option 1: use enum directly)

Change kernel signature to accept the enum:

 inline void convective_u_face_kernel_staggered(
     int i, int j,
     int u_stride, int v_stride, int conv_stride,
-    double dx, double dy, int scheme,
+    double dx, double dy, ConvectiveScheme scheme,
     const double* u_ptr, const double* v_ptr,
     double* conv_u_ptr)
 {
     // ...
-    const bool use_central = (scheme == 0 || scheme == 2);   // 0: Central, 1: Upwind, 2: SkewSymmetric
-    const bool use_skew    = (scheme == 2);
+    const bool use_central = (scheme == ConvectiveScheme::Central || scheme == ConvectiveScheme::SkewSymmetric);
+    const bool use_skew    = (scheme == ConvectiveScheme::SkewSymmetric);

Then at call site:

-    const int scheme = static_cast<int>(config_.convective_scheme);
+    const ConvectiveScheme scheme = config_.convective_scheme;
🔎 Proposed fix (option 2: named constants)

Define constants at file scope:

namespace {
    constexpr int SCHEME_CENTRAL = 0;
    constexpr int SCHEME_UPWIND = 1;
    constexpr int SCHEME_SKEW_SYMMETRIC = 2;
}

Then use in kernel:

-    const bool use_central = (scheme == 0 || scheme == 2);   // 0: Central, 1: Upwind, 2: SkewSymmetric
-    const bool use_skew    = (scheme == 2);
+    const bool use_central = (scheme == SCHEME_CENTRAL || scheme == SCHEME_SKEW_SYMMETRIC);
+    const bool use_skew    = (scheme == SCHEME_SKEW_SYMMETRIC);

Also applies to: 468-469

📜 Review details

Configuration used: defaults

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between bdd6b0b and ee38b5c.

📒 Files selected for processing (2)
  • src/config.cpp
  • src/solver.cpp
🧰 Additional context used
📓 Path-based instructions (2)
**/*.{cpp,cc,cxx}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx}: Do not add platform-specific tolerances with #ifdef __APPLE__ or similar platform checks in tests
Include code comments explaining WHY, not WHAT; document numerical algorithms; note non-obvious optimizations

Files:

  • src/config.cpp
  • src/solver.cpp
**/*.{cpp,cc,cxx,h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx,h,hpp}: Avoid overly strict floating-point comparisons (e.g., == for doubles); use appropriate tolerances based on algorithm
Fix all compiler warnings before pushing; use [[maybe_unused]] for intentionally unused variables in assertions
Use const for read-only references and mark methods const if they don't modify state
Use RAII for resource management; use smart pointers or RAII wrappers instead of manual new/delete or malloc/free
Use OmpDeviceBuffer wrapper for GPU buffer management instead of manual allocation
Check return values and use exceptions for error conditions with informative error messages

Files:

  • src/config.cpp
  • src/solver.cpp
🧠 Learnings (3)
📓 Common learnings
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability

Applied to files:

  • src/solver.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Use `solve_steady_with_snapshots()` for automatic VTK output; specify `num_snapshots` in config; files are numbered sequentially plus final

Applied to files:

  • src/solver.cpp
🧬 Code graph analysis (1)
src/config.cpp (1)
include/config.hpp (1)
  • argc (113-113)
⏰ Context from checks skipped due to timeout of 90000ms. You can increase the timeout in your CodeRabbit configuration to a maximum of 15 minutes (900000ms). (3)
  • GitHub Check: gpu-tests
  • GitHub Check: build-and-test (ubuntu-latest, Debug)
  • GitHub Check: build-and-test (ubuntu-latest, Release)
🔇 Additional comments (5)
src/config.cpp (1)

122-122: LGTM!

The perturbation_amplitude parameter is loaded consistently with other config parameters using the get_double helper.

src/solver.cpp (4)

422-454: Boundary guard correctly prevents out-of-bounds access.

The guard at line 427 (if (i >= 2)) correctly prevents accessing v_ptr at negative indices when the full skew-symmetric stencil cannot be applied near boundaries. The fallback to pure advective form is a reasonable compromise.

Note that this creates a scheme inconsistency: boundary-adjacent cells use advective form while interior cells use skew-symmetric form. For most applications this is acceptable, but it may affect energy conservation properties near boundaries in high-resolution DNS. Consider documenting this behavior in the function's docstring.

Based on learnings, Navier-Stokes solver validation should verify energy conservation and momentum balance, especially with the new skew-symmetric scheme.


503-537: Boundary guard implementation is consistent with u-momentum kernel.

The guard at line 508 (if (j >= 2)) mirrors the approach in the u-momentum kernel and correctly prevents out-of-bounds access. The fallback to advective form for j < 2 is consistent with the u-momentum implementation.

The same consistency considerations apply: boundary-adjacent cells use a different discretization than interior cells. This is acceptable for most applications but may affect conservation properties near walls.


2473-2484: LGTM! Helpful diagnostic output for adaptive time stepping.

The diagnostic output provides useful information for monitoring adaptive time stepping behavior. The implementation correctly computes the current CFL number and only activates when adaptive_dt is enabled, avoiding unnecessary output in fixed-dt simulations.


2516-2532: LGTM! Useful addition of velocity component scalar fields.

Adding u_velocity and v_velocity as separate scalar fields complements the existing vector field and enables better component-wise visualization in ParaView/VisIt. The implementation correctly uses cell-centered interpolation and follows VTK format conventions.

Copy link

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 2

♻️ Duplicate comments (1)
app/main_channel.cpp (1)

154-160: Consider hoisting mode calculations out of inner loops.

The values kx, w_m, A_m, and phase depend only on the mode index m, not on the grid indices i and j. For better performance, especially on large grids, consider precomputing these values in a separate outer loop over modes and storing them in arrays.

Also applies to: 178-185

🧹 Nitpick comments (3)
run_two_build_test.sh (2)

15-15: Replace hardcoded paths with portable alternatives.

The script uses hardcoded absolute paths specific to the user sbryngelson3, reducing portability. Consider using $SLURM_SUBMIT_DIR (the directory from which the job was submitted) or deriving paths relative to the script location.

🔎 Proposed refactor
-cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn
+PROJECT_DIR="${SLURM_SUBMIT_DIR:-$(pwd)}"
+cd "$PROJECT_DIR" || { echo "Failed to cd to project directory"; exit 1; }

Then later:

-cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build
+cd "$PROJECT_DIR/build" || { echo "Failed to cd to build directory"; exit 1; }

Also applies to: 51-51


21-23: Consider cleaning build_cpu directory for test determinism.

Creating build_cpu with mkdir -p preserves any existing content, which could include stale build artifacts or reference data. For fully deterministic tests, consider removing the directory first.

🔎 Proposed refactor
 # Create CPU-only build directory
+rm -rf build_cpu
 mkdir -p build_cpu
 cd build_cpu || { echo "Failed to cd to build_cpu"; exit 1; }
src/config.cpp (1)

130-131: Consider validating enum values for robustness.

The scheme parsing correctly handles multiple aliases, but there's no validation that the assigned enum value is within the valid range. While unlikely to cause issues with the current controlled enum, adding a check like assert(scheme >= 0 && scheme < 3) would catch future enum changes.

Also applies to: 296-297

📜 Review details

Configuration used: defaults

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between ee38b5c and 9625407.

📒 Files selected for processing (5)
  • app/main_channel.cpp
  • run_two_build_test.sh
  • src/config.cpp
  • src/solver.cpp
  • tests/test_solver_cpu_gpu.cpp
🧰 Additional context used
📓 Path-based instructions (3)
**/*.{cpp,cc,cxx}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx}: Do not add platform-specific tolerances with #ifdef __APPLE__ or similar platform checks in tests
Include code comments explaining WHY, not WHAT; document numerical algorithms; note non-obvious optimizations

Files:

  • src/config.cpp
  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • app/main_channel.cpp
**/*.{cpp,cc,cxx,h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx,h,hpp}: Avoid overly strict floating-point comparisons (e.g., == for doubles); use appropriate tolerances based on algorithm
Fix all compiler warnings before pushing; use [[maybe_unused]] for intentionally unused variables in assertions
Use const for read-only references and mark methods const if they don't modify state
Use RAII for resource management; use smart pointers or RAII wrappers instead of manual new/delete or malloc/free
Use OmpDeviceBuffer wrapper for GPU buffer management instead of manual allocation
Check return values and use exceptions for error conditions with informative error messages

Files:

  • src/config.cpp
  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • app/main_channel.cpp
**/*.{sh,slurm,batch}

📄 CodeRabbit inference engine (.cursor/rules/lean-changes.mdc)

Always submit slurm jobs using the embers QOS

Files:

  • run_two_build_test.sh
🧠 Learnings (9)
📓 Common learnings
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to **/*.{cpp,cc,cxx,h,hpp} : Check return values and use exceptions for error conditions with informative error messages

Applied to files:

  • src/config.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • app/main_channel.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_tg_validation.cpp : Taylor-Green Vortex test must validate energy decay with <5% error over 100 timesteps using viscous terms, time integration, and periodic BCs

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: For GPU-related changes, also run `./test_before_ci_gpu.sh` to run the complete GPU CI test suite locally

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
  • run_two_build_test.sh
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Explicitly upload data to GPU; check `USE_GPU_OFFLOAD` is defined; verify `omp_get_num_devices() > 0` at runtime for GPU offload code

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Run `./test_before_ci.sh` before every push to repository to test Debug build, Release build, and all unit tests

Applied to files:

  • run_two_build_test.sh
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Use `solve_steady_with_snapshots()` for automatic VTK output; specify `num_snapshots` in config; files are numbered sequentially plus final

Applied to files:

  • src/solver.cpp
  • app/main_channel.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Always call `set_body_force()` for driven flows; initialize velocity field before solving; set turbulence model before first step

Applied to files:

  • app/main_channel.cpp
🧬 Code graph analysis (3)
src/config.cpp (1)
include/config.hpp (1)
  • argc (113-113)
tests/test_solver_cpu_gpu.cpp (5)
include/turbulence_nn_tbnn.hpp (6)
  • mesh (37-45)
  • mesh (63-63)
  • mesh (106-106)
  • mesh (108-108)
  • mesh (113-120)
  • mesh (123-123)
include/turbulence_transport.hpp (4)
  • mesh (104-111)
  • mesh (139-146)
  • bc (253-253)
  • bc (253-253)
include/solver.hpp (2)
  • bc (62-62)
  • step (193-193)
include/fields.hpp (16)
  • i (19-19)
  • i (19-19)
  • i (20-20)
  • i (20-20)
  • i (62-62)
  • i (62-62)
  • i (63-63)
  • i (63-63)
  • i (67-67)
  • i (67-67)
  • i (68-68)
  • i (68-68)
  • i (86-86)
  • i (95-95)
  • i (96-96)
  • i (116-116)
src/solver.cpp (2)
  • step (1349-1909)
  • step (1349-1349)
app/main_channel.cpp (2)
include/mesh.hpp (4)
  • j (74-74)
  • j (74-74)
  • j (78-78)
  • j (78-78)
include/fields.hpp (12)
  • i (19-19)
  • i (19-19)
  • i (20-20)
  • i (20-20)
  • i (62-62)
  • i (62-62)
  • i (63-63)
  • i (63-63)
  • i (67-67)
  • i (67-67)
  • i (68-68)
  • i (68-68)
🪛 Shellcheck (0.11.0)
run_two_build_test.sh

[warning] 15-15: Use 'cd ... || exit' or 'cd ... || return' in case cd fails.

(SC2164)


[warning] 23-23: Use 'cd ... || exit' or 'cd ... || return' in case cd fails.

(SC2164)


[warning] 51-51: Use 'cd ... || exit' or 'cd ... || return' in case cd fails.

(SC2164)

⏰ Context from checks skipped due to timeout of 90000ms. You can increase the timeout in your CodeRabbit configuration to a maximum of 15 minutes (900000ms). (3)
  • GitHub Check: build-and-test (ubuntu-latest, Release)
  • GitHub Check: build-and-test (ubuntu-latest, Debug)
  • GitHub Check: gpu-tests
🔇 Additional comments (12)
run_two_build_test.sh (1)

1-12: LGTM: SLURM configuration is correct.

The SLURM directives properly use the embers QOS as required by coding guidelines, and the resource allocation (1 GPU, 8 CPUs, 30 minutes) is appropriate for this two-build test workflow.

src/config.cpp (1)

190-214: LGTM! Robust numeric parsing with clear error messages.

The helper functions properly catch both std::invalid_argument and std::out_of_range, providing informative error messages that include the argument name and invalid value before exiting. This addresses the exception-handling concern from previous reviews.

tests/test_solver_cpu_gpu.cpp (3)

288-346: LGTM! Comprehensive skew-symmetric test case.

The new test case properly configures the skew-symmetric scheme, runs the simulation, and collects metrics consistently with the existing test cases. The GPU sync is correctly placed before metrics computation.


485-488: LGTM! Critical GPU syncs properly placed.

The sync_from_gpu() calls are correctly placed after the simulation runs and before comparing CPU and GPU results. The comments clearly explain why these syncs are critical - without them, comparisons would use stale host data instead of the actual GPU results.

Also applies to: 544-547, 619-622, 679-682


558-631: LGTM! Skew-symmetric tests provide good coverage.

The new test functions validate the skew-symmetric scheme for both Taylor-Green vortex (fully periodic) and channel flow (periodic-x, wall-y) configurations. They follow the established test structure and include proper GPU synchronization before comparisons.

Also applies to: 633-691

app/main_channel.cpp (3)

118-193: LGTM! Well-designed multi-mode perturbation field.

The multi-mode streamfunction approach correctly generates a divergence-free velocity field with:

  • 1/k² spectral weighting for realistic turbulence transition
  • Proper amplitude normalization via sqrt(w_m)
  • Wall-compatible factors (sin²) that vanish at y=±1
  • Deterministic phase offsets for mode diversity

The implementation ensures exact incompressibility via the streamfunction formulation (u = ∂ψ/∂y, v = -∂ψ/∂x).


342-342: Critical fix: Adaptive dt now properly applied.

The change from (void)solver.compute_adaptive_dt() to solver.set_dt(solver.compute_adaptive_dt()) is essential - the previous code computed the adaptive time step but never applied it, effectively disabling adaptive time-stepping. This fix ensures the computed dt is actually used by the solver.


311-327: LGTM! Clear base flow initialization logic.

The conditional initialization based on perturbation_amplitude is clean and well-structured:

  • When amplitude > 0: applies multi-mode perturbation on top of base flow
  • Otherwise: starts from uniform base flow
  • Clear logging helps users understand the chosen initialization path
src/solver.cpp (4)

380-380: LGTM! Cleaner scheme parameter design.

The change from a boolean flag to an integer scheme parameter (0=Central, 1=Upwind, 2=SkewSymmetric) is more maintainable and extensible. The local derivation of use_central and use_skew within kernels keeps the logic clear, and all call sites have been consistently updated.

Also applies to: 461-461, 394-395, 469-470, 1046-1046, 1060-1060, 1067-1067, 1075-1075, 1082-1083


422-454: LGTM! Boundary checks prevent out-of-bounds access.

The skew-symmetric kernels now include proper boundary checks:

  • u-momentum: requires i >= 2 to safely access v_ptr[i-2]
  • v-momentum: requires j >= 2 to safely access u_ptr[j-2]

When these conditions aren't met (near boundaries), the code falls back to the standard advective form. This addresses the critical out-of-bounds concern from previous reviews and allows the skew-symmetric scheme to work with the default Nghost=1.

Also applies to: 504-539


2475-2486: LGTM! Useful diagnostics for adaptive time-stepping.

The VTK output now includes diagnostic information showing the current dt and computed CFL number when adaptive time-stepping is enabled. This helps users understand time-step evolution without being overly verbose (only printed when writing VTK files).


2518-2534: LGTM! Additional velocity scalar fields enhance visualization.

Adding separate u_velocity and v_velocity scalar fields complements the existing velocity vector field. This is useful for visualization tools that work better with scalar fields or when analyzing individual velocity components.


module load nvhpc

cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🔴 Critical

Add error handling to cd commands.

All three cd commands lack error handling. If any directory change fails, the script continues executing commands in the wrong directory, potentially corrupting data or producing misleading errors.

🔎 Proposed fix
-cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn
+cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn || { echo "Failed to cd to project directory"; exit 1; }
-cd build_cpu
+cd build_cpu || { echo "Failed to cd to build_cpu"; exit 1; }
-cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build
+cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build || { echo "Failed to cd to build directory"; exit 1; }

Also applies to: 23-23, 51-51

🧰 Tools
🪛 Shellcheck (0.11.0)

[warning] 15-15: Use 'cd ... || exit' or 'cd ... || return' in case cd fails.

(SC2164)

🤖 Prompt for AI Agents
In run_two_build_test.sh around line 15 (and also apply same fix to lines 23 and
51), the bare cd commands lack error handling; update each cd to verify it
succeeded and abort the script if it failed by checking the command's exit
status, printing a clear error message indicating which directory change failed
and exiting with a non‑zero status (or enable strict error handling at the top
with set -e and ensure cd failures will stop execution).

Comment on lines +46 to +52
echo ""
echo "=========================================="
echo "Step 3: Running GPU build to compare"
echo "=========================================="

cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build
./test_solver_cpu_gpu --compare-prefix ../build_cpu/cpu_ref
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟠 Major

Validate GPU build prerequisites before comparison.

The script assumes the GPU build directory exists and contains a compiled test_solver_cpu_gpu binary. If either assumption is false, the test fails with unclear error messages.

🔎 Proposed fix
 echo ""
 echo "=========================================="
 echo "Step 3: Running GPU build to compare"
 echo "=========================================="
 
+GPU_BUILD_DIR="/storage/home/hcoda1/6/sbryngelson3/cfd-nn/build"
+if [ ! -d "$GPU_BUILD_DIR" ]; then
+    echo "GPU build directory does not exist: $GPU_BUILD_DIR"
+    exit 1
+fi
+
+if [ ! -x "$GPU_BUILD_DIR/test_solver_cpu_gpu" ]; then
+    echo "GPU test binary not found or not executable: $GPU_BUILD_DIR/test_solver_cpu_gpu"
+    exit 1
+fi
+
-cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build
+cd "$GPU_BUILD_DIR" || { echo "Failed to cd to GPU build directory"; exit 1; }
 ./test_solver_cpu_gpu --compare-prefix ../build_cpu/cpu_ref
📝 Committable suggestion

‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.

Suggested change
echo ""
echo "=========================================="
echo "Step 3: Running GPU build to compare"
echo "=========================================="
cd /storage/home/hcoda1/6/sbryngelson3/cfd-nn/build
./test_solver_cpu_gpu --compare-prefix ../build_cpu/cpu_ref
echo ""
echo "=========================================="
echo "Step 3: Running GPU build to compare"
echo "=========================================="
GPU_BUILD_DIR="/storage/home/hcoda1/6/sbryngelson3/cfd-nn/build"
if [ ! -d "$GPU_BUILD_DIR" ]; then
echo "GPU build directory does not exist: $GPU_BUILD_DIR"
exit 1
fi
if [ ! -x "$GPU_BUILD_DIR/test_solver_cpu_gpu" ]; then
echo "GPU test binary not found or not executable: $GPU_BUILD_DIR/test_solver_cpu_gpu"
exit 1
fi
cd "$GPU_BUILD_DIR" || { echo "Failed to cd to GPU build directory"; exit 1; }
./test_solver_cpu_gpu --compare-prefix ../build_cpu/cpu_ref
🧰 Tools
🪛 Shellcheck (0.11.0)

[warning] 51-51: Use 'cd ... || exit' or 'cd ... || return' in case cd fails.

(SC2164)

🤖 Prompt for AI Agents
In run_two_build_test.sh around lines 46 to 52, the script assumes the GPU build
directory and the test_solver_cpu_gpu binary exist before running the
comparison; add explicit checks that the directory
/storage/home/hcoda1/6/sbryngelson3/cfd-nn/build exists and is readable, verify
the file ./test_solver_cpu_gpu exists and is executable, and verify the
reference directory ../build_cpu/cpu_ref exists; if any check fails, print a
clear error message identifying which prerequisite is missing and exit with a
non-zero status so the failure is immediate and informative.

Copy link

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 0

Caution

Some comments are outside the diff and can’t be posted inline due to platform limitations.

⚠️ Outside diff range comments (1)
tests/test_solver_cpu_gpu.cpp (1)

693-760: Critical: Missing GPU synchronization before comparison.

The test compares CPU and GPU velocity fields (lines 748-751) but lacks the required GPU-to-host synchronization before the comparison loop. Without solver_gpu.sync_from_gpu(), the host will access stale GPU data, causing incorrect validation results.

🔎 Proposed fix: Add GPU sync after time-stepping loop
         // Run 5 steps
         for (int step = 0; step < 5; ++step) {
             solver_cpu.step();
             solver_gpu.step();
         }
         
+#ifdef USE_GPU_OFFLOAD
+        solver_gpu.sync_from_gpu();
+#endif
+        
         // Quick comparison
         double max_diff = 0.0;
📜 Review details

Configuration used: defaults

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between 9625407 and 2be8c16.

📒 Files selected for processing (1)
  • tests/test_solver_cpu_gpu.cpp
🧰 Additional context used
📓 Path-based instructions (2)
**/*.{cpp,cc,cxx}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx}: Do not add platform-specific tolerances with #ifdef __APPLE__ or similar platform checks in tests
Include code comments explaining WHY, not WHAT; document numerical algorithms; note non-obvious optimizations

Files:

  • tests/test_solver_cpu_gpu.cpp
**/*.{cpp,cc,cxx,h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx,h,hpp}: Avoid overly strict floating-point comparisons (e.g., == for doubles); use appropriate tolerances based on algorithm
Fix all compiler warnings before pushing; use [[maybe_unused]] for intentionally unused variables in assertions
Use const for read-only references and mark methods const if they don't modify state
Use RAII for resource management; use smart pointers or RAII wrappers instead of manual new/delete or malloc/free
Use OmpDeviceBuffer wrapper for GPU buffer management instead of manual allocation
Check return values and use exceptions for error conditions with informative error messages

Files:

  • tests/test_solver_cpu_gpu.cpp
🧠 Learnings (5)
📓 Common learnings
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_tg_validation.cpp : Taylor-Green Vortex test must validate energy decay with <5% error over 100 timesteps using viscous terms, time integration, and periodic BCs

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: For GPU-related changes, also run `./test_before_ci_gpu.sh` to run the complete GPU CI test suite locally

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Explicitly upload data to GPU; check `USE_GPU_OFFLOAD` is defined; verify `omp_get_num_devices() > 0` at runtime for GPU offload code

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
🧬 Code graph analysis (1)
tests/test_solver_cpu_gpu.cpp (5)
include/solver.hpp (1)
  • bc (62-62)
app/main_channel.cpp (1)
  • solver (258-258)
tests/test_time_history_consistency.cpp (2)
  • solver (157-157)
  • solver_gpu (227-227)
tests/test_divergence_all_bcs.cpp (5)
  • solver (66-66)
  • solver (146-146)
  • solver (196-196)
  • solver (245-245)
  • solver (393-393)
app/compare_channel_cpu_gpu.cpp (1)
  • solver_gpu (263-263)
⏰ Context from checks skipped due to timeout of 90000ms. You can increase the timeout in your CodeRabbit configuration to a maximum of 15 minutes (900000ms). (3)
  • GitHub Check: build-and-test (ubuntu-latest, Debug)
  • GitHub Check: build-and-test (ubuntu-latest, Release)
  • GitHub Check: gpu-tests
🔇 Additional comments (5)
tests/test_solver_cpu_gpu.cpp (5)

288-346: LGTM! Well-structured skew-symmetric validation case.

The new Case D properly validates the skew-symmetric convective scheme with Taylor-Green initialization, includes GPU synchronization before metrics collection, and uses a distinct metric prefix for clarity.


485-488: Excellent addition of explicit GPU synchronization.

The GPU sync calls with clear comments ensure host-visible state before CPU/GPU comparisons, preventing subtle test failures from stale data.

Also applies to: 544-547


558-631: LGTM! Comprehensive skew-symmetric Taylor-Green validation.

The test properly exercises the skew-symmetric convective scheme with Taylor-Green vortex initialization, includes GPU synchronization, and validates CPU/GPU consistency.


633-691: LGTM! Channel flow skew-symmetric validation.

The test validates skew-symmetric scheme behavior with channel flow configuration (periodic-x, wall-y boundaries), includes proper GPU sync, and ensures CPU/GPU consistency.


820-826: Clear usage documentation.

The updated usage message properly guides users to the two-build dump/compare workflow for CPU vs GPU validation.

Copy link

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 0

♻️ Duplicate comments (1)
examples/channel/cases/turbulent.cfg (1)

51-51: Increase poisson_max_iter to at least 10,000 for DNS production.

poisson_max_iter = 5 is critically insufficient for a 256×512 DNS turbulent channel run. Inadequate Poisson solver convergence violates the divergence-free constraint and compromises numerical stability. All other production channel/RANS configs use 10,000–20,000 iterations. Use at least 10,000 iterations or justify this as an intentional performance study with clear documentation.

🔎 Recommended fix
-poisson_max_iter = 5
+poisson_max_iter = 10000
🧹 Nitpick comments (2)
examples/channel/cases/turbulent.cfg (2)

17-18: Remove commented duplicate nu definition.

Line 18 duplicates the active line 17 with a different value. If this is exploration of different Reynolds numbers, document the choice in a comment rather than leaving dead code.

🔎 Cleanup suggestion
 nu = 0.00666666   # = 1/300 (higher Re_tau)
-# nu = 0.0033333333333333333   # = 1/300 (higher Re_tau)

38-39: Remove commented duplicate convective_scheme line.

Line 39 duplicates the active line 38. Clean up experimental artifacts to improve config readability.

🔎 Cleanup suggestion
 convective_scheme = skew_symmetric # Energy-conserving for DNS
-# convective_scheme = skew_symmetric  # Energy-conserving for DNS
📜 Review details

Configuration used: defaults

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between 2be8c16 and 80371c5.

📒 Files selected for processing (3)
  • examples/channel/cases/turbulent.cfg
  • src/solver.cpp
  • tests/test_solver_cpu_gpu.cpp
🧰 Additional context used
📓 Path-based instructions (2)
**/*.{cpp,cc,cxx}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx}: Do not add platform-specific tolerances with #ifdef __APPLE__ or similar platform checks in tests
Include code comments explaining WHY, not WHAT; document numerical algorithms; note non-obvious optimizations

Files:

  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
**/*.{cpp,cc,cxx,h,hpp}

📄 CodeRabbit inference engine (.cursorrules)

**/*.{cpp,cc,cxx,h,hpp}: Avoid overly strict floating-point comparisons (e.g., == for doubles); use appropriate tolerances based on algorithm
Fix all compiler warnings before pushing; use [[maybe_unused]] for intentionally unused variables in assertions
Use const for read-only references and mark methods const if they don't modify state
Use RAII for resource management; use smart pointers or RAII wrappers instead of manual new/delete or malloc/free
Use OmpDeviceBuffer wrapper for GPU buffer management instead of manual allocation
Check return values and use exceptions for error conditions with informative error messages

Files:

  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
🧠 Learnings (10)
📓 Common learnings
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: For GPU-related changes, also run `./test_before_ci_gpu.sh` to run the complete GPU CI test suite locally

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_physics_validation.cpp : Navier-Stokes solver validation tests must verify: viscous terms, pressure gradient, parabolic velocity profiles, divergence-free constraint, momentum balance, channel symmetry, cross-model consistency, and numerical stability

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
  • src/solver.cpp
  • examples/channel/cases/turbulent.cfg
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to **/*.{cpp,cc,cxx,h,hpp} : Use `OmpDeviceBuffer` wrapper for GPU buffer management instead of manual allocation

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to **/*.{cpp,cc,cxx} : Do not add platform-specific tolerances with `#ifdef __APPLE__` or similar platform checks in tests

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Applies to tests/test_tg_validation.cpp : Taylor-Green Vortex test must validate energy decay with <5% error over 100 timesteps using viscous terms, time integration, and periodic BCs

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Tests should pass in both Debug and Release builds without modification

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Explicitly upload data to GPU; check `USE_GPU_OFFLOAD` is defined; verify `omp_get_num_devices() > 0` at runtime for GPU offload code

Applied to files:

  • tests/test_solver_cpu_gpu.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Use `solve_steady_with_snapshots()` for automatic VTK output; specify `num_snapshots` in config; files are numbered sequentially plus final

Applied to files:

  • src/solver.cpp
📚 Learning: 2025-12-25T20:38:43.379Z
Learnt from: CR
Repo: comp-physics/cfd-nn PR: 0
File: .cursorrules:0-0
Timestamp: 2025-12-25T20:38:43.379Z
Learning: Always call `set_body_force()` for driven flows; initialize velocity field before solving; set turbulence model before first step

Applied to files:

  • examples/channel/cases/turbulent.cfg
🧬 Code graph analysis (1)
src/solver.cpp (1)
include/solver.hpp (3)
  • dp_dx (212-212)
  • pressure_ (106-106)
  • pressure_ (113-113)
⏰ Context from checks skipped due to timeout of 90000ms. You can increase the timeout in your CodeRabbit configuration to a maximum of 15 minutes (900000ms). (3)
  • GitHub Check: gpu-tests
  • GitHub Check: build-and-test (ubuntu-latest, Debug)
  • GitHub Check: build-and-test (ubuntu-latest, Release)
🔇 Additional comments (7)
tests/test_solver_cpu_gpu.cpp (2)

288-346: LGTM: Skew-symmetric scheme test case.

The test properly validates the new SkewSymmetric convective scheme using a Taylor-Green vortex, consistent with the existing test structure. Metrics collection and GPU synchronization are correctly implemented.


410-417: LGTM: Clear two-build workflow documentation.

The usage message accurately describes the new dump/compare workflow and explicitly states that legacy single-binary tests have been removed. This will help users understand the testing approach.

src/solver.cpp (5)

458-540: LGTM: v-momentum skew-symmetric kernel with proper boundary checks.

The kernel correctly mirrors the u-momentum implementation with appropriate boundary checks (j >= 2 and (j+1) <= u_max_j) to prevent out-of-bounds access when computing kinetic energy gradients. The fallback to pure advective form near boundaries maintains correctness.


1212-1265: LGTM: Pressure gradient computation with periodic boundary handling.

The implementation correctly handles periodic boundaries in the x-direction by wrapping to the opposite edge when computing dp/dx at boundary cells. The y-direction uses one-sided differences at boundaries, appropriate for channel flow with no-slip walls. The logic aligns with the staggered grid discretization.


2518-2529: LGTM: VTK diagnostic output for adaptive time stepping.

The diagnostic correctly reports the current dt and computed CFL when writing VTK files during adaptive time-stepping runs. This is valuable for monitoring the solver's time-step selection and verifying CFL constraint satisfaction in DNS simulations.


2561-2601: LGTM: Extended VTK output with velocity components and pressure gradient.

The additions provide valuable diagnostic fields for post-processing:

  • Scalar fields for u and v velocity components (interpolated to cell centers)
  • Pressure gradient vector field computed via central differences

The const_cast on line 2593 is necessary because compute_pressure_gradient modifies temporary scratch fields dp_dx and dp_dy, even though it doesn't modify persistent solver state. This is acceptable for VTK output.


377-455: Boundary check for skew-symmetric kernel is correct—v_stride properly validates column indices.

The check (i + 1) < v_stride correctly uses v_stride = Nx + 2*Ng as the upper bound for v-velocity column access. Combined with i >= 2, these conditions prevent out-of-bounds access when reading v-velocities at neighboring i±1 positions. The fallback to pure advective form near boundaries ensures correctness. No changes needed.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Development

Successfully merging this pull request may close these issues.

2 participants