Module simd_ops 

SIMD-accelerated operations for QSM processing

This module provides vectorized versions of common operations used in iterative algorithms like MEDI, TV, and TGV. When the simd feature is enabled, these use 128-bit SIMD (f32x4) which is compatible with both native SSE/NEON and WASM SIMD.

All operations have scalar fallbacks when SIMD is disabled.

Constants

SIMD_WIDTH

SIMD lane width (4 for f32x4)

Functions

apply_gradient_weights_f32

Apply per-direction gradient weights: out[i] = mx[i] * p[i] * mx[i] * gx[i] This is the core operation in MEDI’s regularization term

apply_gradient_weights_f64

Apply per-direction gradient weights (f64 version)

axpy_f32

Compute a[i] = a[i] + alpha * b[i] (axpy operation)

axpy_f64

Compute a[i] = a[i] + alpha * b[i] (axpy operation, f64 version)

combine_terms_f32

Combine regularization and data terms: out[i] = lambda * reg[i] + data[i] Matches MATLAB MEDI: y = D + R where D is data term, R = lambda * reg term

compute_p_weights_f32

Compute P = 1 / sqrt(ux^2 + uy^2 + uz^2 + beta) where ux = mx * gx, uy = my * gy, uz = mz * gz

dot_product_f32

Compute dot product: sum(a[i] * b[i])

dot_product_f64

Compute dot product: sum(a[i] * b[i]) (f64 version)

negate_f32

Negate array in place: a[i] = -a[i]

negate_f64

Negate array in place: a[i] = -a[i] (f64 version)

norm_squared_f32

Compute squared norm: sum(a[i]^2)

norm_squared_f64

Compute squared norm: sum(a[i]^2) (f64 version)

scale_f32

Scale array in place: a[i] = alpha * a[i]

scale_f64

Scale array in place: a[i] = alpha * a[i] (f64 version)

subtract_f32

Subtract arrays element-wise: out[i] = a[i] - b[i]

subtract_f64

Subtract arrays element-wise: out[i] = a[i] - b[i] (f64 version)

xpby_f32

Compute a[i] = b[i] + beta * a[i] (used in CG for p update)