Module l2_norm 

L2 norm using BLAS/LAPACK-style incremental scaling.

L2 Norm (Euclidean Norm) Computation

This module provides a numerically stable implementation of the L2 norm (Euclidean norm) for slices of RealScalar values.

Algorithm: BLAS/LAPACK-Style Incremental Scaling

The implementation uses the scaled sum of squares approach, which is the same algorithm used in modern BLAS/LAPACK implementations (e.g., DNRM2). This single-pass algorithm avoids overflow and underflow by maintaining the invariant:

$$|x|_2 = \text{scale} \cdot \sqrt{\text{sumsq}}$$

where:

• scale is the maximum absolute value seen so far
• sumsq is the sum of squared normalized values: $\sum_i (|x_i| / \text{scale})^2$

Key Properties

Property	Description
Single-pass	O(n) time complexity with one iteration over the data
No overflow	All squared values are ≤ 1 (normalized by max)
No underflow	Scale factor preserves the magnitude of small values
Order-independent	Same result regardless of element ordering
Backend-agnostic	Works with both f64 and rug::Float backends

How It Works

1. Track scale (max |xᵢ| seen) and sumsq (sum of normalized squares)
2. For each element xᵢ:
   • If |xᵢ| > scale: rescale sumsq to new scale, then accumulate
   • Otherwise: just accumulate (|xᵢ|/scale)²
3. Return scale × √sumsq

Why Not Naive Implementation?

The naive approach sqrt(sum(x²)) fails for extreme values:

Naive overflow:   ||[1e200, 1e200]|| → (1e200)² = Inf → sqrt(Inf) = Inf  ❌
Naive underflow:  ||[1e-200, 1e-200]|| → (1e-200)² = 0 → sqrt(0) = 0    ❌

Scaled approach:  ||[1e200, 1e200]|| → 1e200 × sqrt(2) ≈ 1.41e200        ✅
Scaled approach:  ||[1e-200, 1e-200]|| → 1e-200 × sqrt(2) ≈ 1.41e-200    ✅

Usage Examples

Basic Usage

use num_valid::{RealNative64StrictFinite, RealScalar, algorithms::l2_norm::l2_norm};

let data: Vec<RealNative64StrictFinite> = vec![
    RealNative64StrictFinite::from_f64(3.0),
    RealNative64StrictFinite::from_f64(4.0),
];

let norm = l2_norm(&data);
assert_eq!(*norm.as_ref(), 5.0); // 3² + 4² = 25, √25 = 5

Handling Extreme Values

use num_valid::{RealNative64StrictFinite, RealScalar, algorithms::l2_norm::l2_norm};

// Values near overflow threshold - naive approach would fail
let large: Vec<RealNative64StrictFinite> = vec![
    RealNative64StrictFinite::from_f64(1e154),
    RealNative64StrictFinite::from_f64(1e154),
];
let norm = l2_norm(&large);
// Result: √2 × 1e154 ≈ 1.41e154 (no overflow!)
assert!((norm.as_ref() / 1e154 - std::f64::consts::SQRT_2).abs() < 1e-10);

// Values near underflow threshold
let small: Vec<RealNative64StrictFinite> = vec![
    RealNative64StrictFinite::from_f64(1e-154),
    RealNative64StrictFinite::from_f64(1e-154),
];
let norm = l2_norm(&small);
// Result: √2 × 1e-154 ≈ 1.41e-154 (no underflow!)
assert!((norm.as_ref() / 1e-154 - std::f64::consts::SQRT_2).abs() < 1e-10);

With Arbitrary-Precision Backend

use num_valid::{RealRugStrictFinite, algorithms::l2_norm::l2_norm};
use try_create::TryNew;

type R = RealRugStrictFinite<200>; // 200-bit precision

// Values beyond f64 range (would be infinity in f64)
let huge = R::try_new(rug::Float::with_val(200, rug::Float::parse("1e1000")?))?;
let data: Vec<R> = vec![huge.clone(), huge.clone()];

let norm = l2_norm(&data);
// Result: √2 × 1e1000 - computed exactly with arbitrary precision

Edge Cases

Input	Result
Empty slice []	0
Single element [x]	`
All zeros [0, 0, 0]	0
Contains zeros [0, 3, 0, 4]	Same as [3, 4] → 5
Negative values [-3, -4]	Same as [3, 4] → 5

Performance Characteristics

• Time complexity: O(n) - single pass over data
• Space complexity: O(1) - only two accumulator variables
• Operations per element: 1 comparison, 1-2 divisions, 1 multiply-add

The algorithm performs more divisions than the naive approach, but this tradeoff is necessary for numerical stability. For performance-critical code where values are known to be in a safe range, consider the naive approach with appropriate bounds checking.

References

• Blue, J. L. (1978). “A Portable Fortran Program to Find the Euclidean Norm of a Vector”. ACM Transactions on Mathematical Software, 4(1), 15-23.
• Anderson, E. (2017). “Algorithm 978: Safe Scaling in the Level 1 BLAS”. ACM Transactions on Mathematical Software, 44(1), 1-28.
• LAPACK Working Note 148: “On Computing LAPACK’s XNRM2”

Functions

l2_norm

Computes the L2 norm (Euclidean norm) of a slice of real scalars.