Constant BOUNDARY_ENFORCE_WGSL 

pub const BOUNDARY_ENFORCE_WGSL: &str = r#"
pub(super) struct BoundaryParams {
    n: u32,
    box_min_x: f32,
    box_min_y: f32,
    box_min_z: f32,
    box_max_x: f32,
    box_max_y: f32,
    box_max_z: f32,
    restitution: f32,
}

@group(0) @binding(0) var<storage, read_write> positions:  array<vec4<f32>>;
@group(0) @binding(1) var<storage, read_write> velocities: array<vec4<f32>>;
@group(0) @binding(2) var<uniform>             params:     BoundaryParams;

@compute @workgroup_size(64)
pub(super) fn main(@builtin(global_invocation_id) id: vec3<u32>) {
    let i = id.x;
    if (i >= params.n) { return; }

    var pos = positions[i].xyz;
    var vel = velocities[i].xyz;
    let e = params.restitution;

    if (pos.x < params.box_min_x) { pos.x = params.box_min_x; vel.x = -vel.x * e; }
    if (pos.x > params.box_max_x) { pos.x = params.box_max_x; vel.x = -vel.x * e; }
    if (pos.y < params.box_min_y) { pos.y = params.box_min_y; vel.y = -vel.y * e; }
    if (pos.y > params.box_max_y) { pos.y = params.box_max_y; vel.y = -vel.y * e; }
    if (pos.z < params.box_min_z) { pos.z = params.box_min_z; vel.z = -vel.z * e; }
    if (pos.z > params.box_max_z) { pos.z = params.box_max_z; vel.z = -vel.z * e; }

    positions[i]  = vec4<f32>(pos, positions[i].w);
    velocities[i] = vec4<f32>(vel, 0.0);
}
"#;

WGSL compute shader for boundary condition enforcement.