OTOM SPH

High-Performance Fluid Dynamics & Volumetric Visualization

Technical Whitepaper v2.0

A Professional Reference for Engineering Validation and Industry-Grade Simulations

OTOMcomposite Engineering Suite
www.otomcomposite.eu | info@otomcomposite.eu

1. Executive Summary

The OTOM SPH program is a cutting-edge Computational Fluid Dynamics (CFD) platform designed for real-time visualization and preliminary engineering analysis. Unlike traditional desktop-bound CFD suites, OTOM SPH leverages the WebGPU API to bring workstation-class simulation performance directly to the web browser.

Key Innovation: By combining Weakly Compressible SPH (WCSPH) with Sub-Particle Scale (SPS) turbulence modeling, the engine provides mathematically grounded results suitable for aerodynamics validation, chemical mixing analysis, and bio-reactor design.

This document details the mathematical framework, the GPU architecture optimizations (including fixed-point atomic operations), and the four distinct rendering pipelines that allow for both raw data visualization and photo-realistic fluid rendering.

2. Physical Foundations of SPH

Smoothed Particle Hydrodynamics (SPH) is a Lagrangian method, meaning the fluid is represented as a collection of discrete moving elements (particles) rather than a fixed grid. This makes it exceptionally well-suited for simulations involving free surfaces, high-velocity impacts, and complex moving boundaries.

2.1 The Navier-Stokes Equations

The simulation solves the conservation of mass and momentum in a Lagrangian frame:

$$ \frac{d\mathbf{v}}{dt} = -\frac{1}{\rho} \nabla p + \nu \nabla^2 \mathbf{v} + \mathbf{g} $$

Where $\rho$ is density, $p$ is pressure, $\nu$ is kinematic viscosity, and $\mathbf{g}$ is gravity.

2.2 The Smoothing Kernel

A continuous field $A(\mathbf{x})$ is approximated by summing the contributions of all particles within a smoothing radius $h$. The engine uses the **Spiky Kernel** for pressure forces to prevent particle clustering:

$$ W_{spiky}(r, h) = \frac{15}{\pi h^6} (h - r)^3 $$

2.3 Kernel Consistency & Normalization

A critical requirement for SPH accuracy is the Normalization Condition. The smoothing kernel $W$ must integrate to unity over its support domain:

$$ \int_{\Omega} W(\mathbf{x} - \mathbf{x}', h) d\mathbf{x}' = 1 $$

In practice, near boundaries or areas with low particle density, this condition is violated (the "Kernel Deficiency" problem). The OTOM SPH engine applies Shepard Renormalization to correct the density field:

$$ \rho_i^{corrected} = \frac{\sum_j m_j W_{ij}}{\sum_j \frac{m_j}{\rho_j} W_{ij}} $$

This ensures that even at the fluid's free surface, the pressure calculations remain stable and do not suffer from unphysical "surface tension" artifacts inherent in basic SPH implementations.

2.4 Time Integration: The Symplectic Euler Scheme

To maintain energy conservation over long simulation periods, the engine uses a Symplectic Euler integrator. Unlike standard Euler, it updates velocity first and then uses the new velocity to update position:

$$ \mathbf{v}_{t+\Delta t} = \mathbf{v}_t + \mathbf{a}_t \Delta t $$ $$ \mathbf{x}_{t+\Delta t} = \mathbf{x}_t + \mathbf{v}_{t+\Delta t} \Delta t $$

This semi-implicit approach significantly improves the stability of the simulation, allowing for larger time-steps ($\Delta t$) without numerical explosion.

3. Advanced Aerodynamic Solver Suite

Industrial Excellence

To move beyond simple "visual water," the OTOM SPH engine includes a suite of advanced solvers that model complex fluid phenomena typically found in high-end engineering software.

3.1 SPS-LES Turbulence Modeling

The Large Eddy Simulation (LES) approach handles sub-grid turbulence. In our SPH context, it models the energy dissipation of vortices smaller than the particle size. This is critical for capturing the "boiling" wake behind an aircraft wing or a rotating blade.

The solver computes the **Rate-of-Strain Tensor** ($S_{ij}$) for every particle:

$$ S_{ij} = \frac{1}{2} \left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) $$

The **Eddy Viscosity** is then injected back into the simulation, preventing the "rubbery" look of basic SPH and providing a physically accurate representation of energy loss in turbulent flows.

3.2 Logarithmic Law of the Wall

Standard simulations often fail at the boundary layer (where the fluid touches an object). We implement the Log-Law of the Wall, which mathematically describes the velocity profile near a surface. This allows the simulation to accurately predict skin friction and drag, even with a relatively low particle count.

4. Real-time CFD Metrics & Analysis

OTOM SPH doesn't just show fluid; it analyzes it. The system provides real-time feedback on critical engineering coefficients, processed directly on the GPU using atomic reduction algorithms.

Metric Engineering Significance Implementation
Drag Coeff (Cd) Measures aerodynamic efficiency. Momentum transfer integral.
Lift Coeff (Cl) Crucial for aerospace/wings. Integrated surface pressure.
Max Shear Rate Cell damage in Bio-Reactors. Fixed-point atomic max on GPU.
Reynolds Number Defines Laminar vs Turbulent. Computed from smoothing radius.
Mixing Efficiency Stirred tank homogeneity. Standard deviation of particle homogeneity.

4.1 High-Precision Force Extraction

The solver calculates the force exerted by the fluid on the obstacle in two concurrent ways:

The Impulse Method: Tracks the exact momentum change ($\Delta \mathbf{p} = m \Delta \mathbf{v}$) of particles as they collide with the surface. For a set of $N_c$ colliding particles: $$ \mathbf{F}_{impulse} = \frac{1}{\Delta t} \sum_{i=1}^{N_c} m_i (\mathbf{v}_{i, post} - \mathbf{v}_{i, pre}) $$ This is highly robust for high-speed flows and "splashy" interactions.
The Pressure Method: Performs a surface integral of the fluid pressure field across the obstacle's geometry. This is more accurate for slow, pressurized flows: $$ \mathbf{F}_{pressure} = \int_{S} p(\mathbf{x}) \mathbf{n} dS \approx \sum_{i} p_i \mathbf{n}_i A_i $$ Where $A_i$ is the effective area of the particle interaction on the surface.

4.2 Engineering Metrics: Detailed Analysis

The OTOM SPH engine extracts high-fidelity data from the particle field. Below is a detailed breakdown of the primary engineering outputs, including their mathematical derivation and practical application examples.

A. Power Consumption (W)

Definition: In the context of Bio-Reactors and mixing tanks, this measures the mechanical power required to drive an impeller through the fluid.

How it's calculated: The engine monitors the rotational torque ($\tau$) acting on the spinning obstacle. Power is derived as $P = \tau \cdot \omega$, where $\omega$ is the angular velocity. The torque is computed by integrating the pressure and viscosity forces across the impeller's surface area.

Real-world Example: An engineer designing a vaccine production tank uses this to size the electric motor. If the simulation shows a power requirement of 450W for a specific impeller design, they might select a 750W motor to ensure a safety margin.

B. Mixing Efficiency (%)

Definition: Quantifies how effectively the fluid is being homogenized. A value of 100% indicates a perfectly uniform mixture.

How it's calculated: The solver tracks a "scalar field" (represented by particle color/age). Efficiency is derived from the Coefficient of Variation (CoV). As the standard deviation of the particle distribution decreases relative to the mean, the efficiency percentage increases.

Real-world Example: In a chemical plant mixing two reactants, 95% efficiency ensures that the reaction occurs uniformly, preventing "hot spots" that could damage the product or the vessel.

C. Wave Resistance (N)

Definition: The component of the fluid's drag force acting on a partially submerged object (like a ship hull) specifically due to the energy spent creating surface waves.

How it's calculated: By integrating the horizontal pressure force acting on the hull's wetted surface area. It is heavily influenced by the Froude Number ($Fr = v / \sqrt{gL}$).

Real-world Example: A naval architect tests a new hull shape. If the Wave Resistance drops from 5000N to 4200N in the simulation, it translates directly to thousands of liters of fuel saved per year for a full-sized vessel.

D. Pressure Drop ($\Delta P$) (Pa)

Definition: The difference in total pressure between the inlet and outlet of a pipe or valve.

How it's calculated: The engine measures the average particle pressure at the "Entry Face" and subtracts the average pressure at the "Exit Face." This takes into account internal friction (viscosity) and turbulence losses.

Real-world Example: In an HVAC system, a high pressure drop across a filter means the fans have to work harder. Engineers use this metric to optimize the internal baffles of a valve to minimize energy loss.

E. Peak Impact Pressure (Pa)

Definition: The maximum instantaneous pressure recorded on a solid boundary during a high-velocity impact event (e.g., a wave hitting a sea wall).

How it's calculated: The solver identifies the particle-to-wall collision with the highest normal force during a specific time-step and divides it by the effective particle area.

Real-world Example: Civil engineers designing a flood barrier for a city use this to determine the thickness of reinforced concrete. If the "Flood-Impact" case shows 2.5 MPa of pressure, the barrier must be built to withstand at least 5.0 MPa for safety.

F. Total Impulse (N·s)

Definition: The integral of the force over time. It represents the total change in momentum transferred from the fluid to an object during an event.

How it's calculated: $J = \int F dt$. The engine accumulates the total force acting on the obstacle over the entire simulation run.

Real-world Example: Used to assess the stability of offshore oil rigs. A single large wave might have a high peak pressure, but the Total Impulse determines if the entire rig will actually shift or tip on its foundation.

5. WebGPU Engine Architecture

The engine is built on a massive-parallel architecture designed to saturate high-end GPU hardware. It utilizes Compute Shaders for 100% of the physics calculations, keeping the CPU free for UI and high-level logic.

5.1 Neighbor Search via Spatial Hashing

For $N$ particles, a naive neighbor search is $O(N^2)$. OTOM SPH uses a Uniform Grid Spatial Hashing algorithm:

$$ \mathbf{c} = \lfloor \mathbf{p} / h \rfloor $$ $$ \mathcal{H}(\mathbf{c}) = (c_x \cdot P_1 \oplus c_y \cdot P_2 \oplus c_z \cdot P_3) \pmod{N_{table}} $$

Where $\mathbf{c}$ is the cell coordinate, $P_{1,2,3}$ are large primes, and $\oplus$ denotes the bitwise XOR operation. By sorting particles by this hash using a GPU Bitonic Sort, we achieve Cache Coherency—neighboring particles in 3D space are usually adjacent in GPU memory, drastically reducing memory latency during the density and viscosity passes.

5.2 Fixed-Point Atomic Operations

GPUs natively struggle with atomic additions on floating-point numbers. OTOM SPH solves this by converting physical units into large integers (fixed-point):

$$ V_{int} = \text{round}(V_{float} \cdot 1000.0) $$ $$ \text{atomicMax}(\text{buffer}[i], V_{int}) $$

This allows for deterministic, race-condition-free global value reduction (like finding the Maximum Shear Rate) across thousands of parallel threads on the GPU hardware.

6. High-Fidelity Rendering Modalities

The system offers four distinct rendering pipelines to suit different needs.

6.1 Screen Space Fluid Rendering (SSFR)

This is the most realistic mode for liquids. It uses a Separable Bilateral Depth Filter to reconstruct a continuous surface from discrete particles. The filter weights $W$ consider both spatial distance and depth similarity to preserve sharp edges:

$$ W(x, y) = G_{\sigma_s}(||x - y||) \cdot G_{\sigma_r}(|D(x) - D(y)|) $$

6.2 Volumetric Raymarching

Ideally suited for gas, smoke, or dense chemicals. The GPU casts rays through a 3D density volume, integrating light transmittance at every step. It features **Physical Absorption (Beer-Lambert)** where the fluid color is a function of the extinction coefficients $\beta_{rgb}$ and the path length $L$:

$$ I(L) = I_0 \cdot e^{-\beta L} $$

7. Industrial Use Cases

Aerospace & Automotive (Wind Tunnel): Rapidly prototype aerodynamic shapes (STL import) and visualize vortex shedding in real-time. The engine calculates Lift and Drag coefficients ($C_l, C_d$) based on surface pressure integration.
Bio-Reactor Engineering: Optimize mixing efficiency in stirred tanks. This use case specifically monitors True Maximum Shear Rate to prevent delicate cell damage and calculates power consumption based on motor torque requirements.
Marine Hydrodynamics (Ship Hull): Analyze flow around partially submerged structures. The system extracts Wave Resistance Force and approximates Wetted Surface Area, critical for naval architecture and fuel efficiency studies.
Process Engineering (Valve/Pipe Flow): Model internal flows through complex plumbing or valves. It calculates the Pressure Drop ($\Delta P$) across obstacles and monitors the Volumetric Flow Rate ($Q$), assisting in pump sizing and baffle optimization.
Civil Engineering (Flood/Impact): Simulate catastrophic events such as dam breaks or urban flooding. The solver extracts Peak Impact Pressure and Total Impulse on infrastructure, enabling safety validation against extreme hydraulic loads.

8. Validation & Accuracy

To ensure engineering relevance, we validated the engine against classic CFD benchmarks:

9. User Operations Manual

The Graphical User Interface (GUI) provides complete, real-time control over the simulation. Below is a detailed breakdown of every menu tab and its associated inputs:

9.1 🌐 Global Settings

9.2 〰️ Streamlines

9.3 ⚙️ Solver

9.4 🧪 Fluid

9.5 🧱 Obstacle

9.6 📥 Container

9.7 ☀️ Environment

9.8 📊 Drag and Lift Coeff

9.9 💼 Use Cases

9.10 🖱️ Interaction

Performace ⚡

11. Detailed Algorithm Implementation

This section provides a low-level look at the WebGPU implementation details that power the OTOM SPH engine.

11.1 The Spiky Smoothing Kernel (WGSL)

The following WGSL snippet demonstrates how the engine efficiently calculates pressure gradients while maintaining numerical stability near the kernel radius boundary:

fn spikyGradient(offset: vec3<f32>, dst: f32, radius: f32, scale: f32) -> vec3<f32> {
    if (dst < radius && dst > 1e-6) {
        let v = radius - dst;
        // Direction points strictly along the radial vector
        return -scale * v * v * (offset / dst);
    }
    return vec3<f32>(0.0);
}
        

11.2 GPU Bitonic Merge Sort

To maintain high cache coherency, particles are sorted by their spatial hash every frame. The engine utilizes a **Bitonic Merge Sort**, a parallel sorting algorithm with $O(N \log^2 N)$ complexity. Unlike CPU-based quicksort, bitonic sort has a fixed execution pattern, making it ideal for the SIMD architecture of modern GPUs.

12. Engineering Setup & Calibration

For reproducible results, the simulation must be calibrated correctly for the specific hardware and fluid type.

12.1 Recommended Hardware

12.2 Stability Calibration

If the simulation becomes unstable (particles "exploding"), follow these engineering steps:

  1. Decrease Time Scale: Reducing the global $dt$ improves solver convergence.
  2. Increase Smoothing Radius: Ensuring each particle has at least 30-40 neighbors prevents pressure spikes.
  3. Tune Pressure Multiplier: In WCSPH, a multiplier that is too high causes the system to become "stiff" and unstable.

10. Glossary & References

WCSPH: Weakly Compressible Smoothed Particle Hydrodynamics.

CFL Condition: Courant-Friedrichs-Lewy condition, ensuring the simulation time-step is small enough for stability.

IOR: Index of Refraction, determining how light bends as it enters the fluid.

References:
1. Monaghan, J.J. (1992). Smoothed Particle Hydrodynamics.
2. Müller, M. (2003). Particle-Based Fluid Simulation for Interactive Applications.
3. Smagorinsky, J. (1963). General Circulation Experiments with the Primitive Equations.

13. Industrial Benchmarks & Performance

The engine is designed for high-throughput engineering validation. Performance scales linearly with the number of GPU cores.

Particle Count Hardware Solver Time (ms) Render Time (ms)
100,000 NVIDIA RTX 3060 1.2 ms 0.8 ms
500,000 NVIDIA RTX 3060 4.5 ms 2.1 ms
1,000,000 NVIDIA RTX 4090 3.8 ms 1.2 ms

13.1 Memory Throughput

By utilizing Unified Memory Architecture in WebGPU, the simulation achieves a theoretical peak bandwidth of ~400 GB/s on high-end desktop GPUs. This allows for near-instantaneous visualization of complex aerodynamic wakes without the latency typically associated with transferring data from the CPU to the display.

© 2026 OTOMcomposite. All mathematical models and software implementations are proprietary.