From 476fecc6b5a8a6a6e800184f57227452fbb8ffbb Mon Sep 17 00:00:00 2001
From: Julian Litz <91479202+julianlitz@users.noreply.github.com>
Date: Fri, 27 Feb 2026 12:55:16 +0100
Subject: [PATCH 1/2] Remove l2_norm_squared and improve comments

Removed the l2_norm_squared function and updated comments for clarity.
---
 .../LinearAlgebra/Vector/vector_operations.h  | 21 +++++--------------
 1 file changed, 5 insertions(+), 16 deletions(-)

diff --git a/include/LinearAlgebra/Vector/vector_operations.h b/include/LinearAlgebra/Vector/vector_operations.h
index 2546a6e0..cd877c58 100644
--- a/include/LinearAlgebra/Vector/vector_operations.h
+++ b/include/LinearAlgebra/Vector/vector_operations.h
@@ -112,18 +112,7 @@ T l1_norm(ConstVector<T> x)
     return result;
 }
 
-template <typename T>
-T l2_norm_squared(ConstVector<T> x)
-{
-    T result      = 0.0;
-    std::size_t n = x.size();
-#pragma omp parallel for reduction(+ : result) if (n > 10'000)
-    for (std::size_t i = 0; i < n; ++i) {
-        result += x(i) * x(i);
-    }
-    return result;
-}
-
+// Underflow- and overflow-resistant implementation of the L2 norm
 template <typename T>
 T l2_norm(ConstVector<T> x)
 {
@@ -137,17 +126,17 @@ T l2_norm(ConstVector<T> x)
             scale = abs_val;
         }
     }
-    if (equals(scale, T{0})) {
+    // 2) if the largest absolute value is zero, the norm is zero
+    if (scale == T{0})
         return T{0};
-    }
-    // 2) accumulate sum of squares of scaled entries
+    // 3) accumulate sum of squares of scaled entries
     T sum = 0.0;
 #pragma omp parallel for reduction(+ : sum) if (n > 10'000)
     for (std::size_t i = 0; i < n; ++i) {
         T value = x(i) / scale;
         sum += value * value;
     }
-    // 3) rescale
+    // 4) rescale
     return scale * std::sqrt(sum);
 }
 

From 5d93c332756a2c9430e34fd8294519c8f9868aa5 Mon Sep 17 00:00:00 2001
From: Julian Litz <91479202+julianlitz@users.noreply.github.com>
Date: Fri, 27 Feb 2026 12:56:08 +0100
Subject: [PATCH 2/2] Refactor l2 norm tests and add overflow/underflow cases

---
 .../Vector/vector_operations.cpp              | 36 ++++++++++++++-----
 1 file changed, 27 insertions(+), 9 deletions(-)

diff --git a/tests/LinearAlgebra/Vector/vector_operations.cpp b/tests/LinearAlgebra/Vector/vector_operations.cpp
index d411a48f..bf63ddb2 100644
--- a/tests/LinearAlgebra/Vector/vector_operations.cpp
+++ b/tests/LinearAlgebra/Vector/vector_operations.cpp
@@ -163,28 +163,46 @@ TEST(VectorOperations, l1_vector_norm)
     EXPECT_DOUBLE_EQ(l1_norm(ConstVector<double>(v)), 8.0);
 }
 
-/* T l2_norm_squared(ConstVector<T>& x); */
+/* T l2_norm(ConstVector<T>& x); */
 
-TEST(VectorOperations, l2_vector_norm_squared)
+TEST(VectorOperations, l2_vector_norm)
 {
     Vector<double> v("v", 3);
     v(0) = 1;
     v(1) = -5;
     v(2) = 2;
     ConstVector<double> const_v(v);
-    EXPECT_DOUBLE_EQ(l2_norm_squared(const_v), 30.0);
+    EXPECT_DOUBLE_EQ(l2_norm(ConstVector<double>(const_v)), std::sqrt(30.0));
 }
 
-/* T l2_norm(ConstVector<T>& x); */
+TEST(VectorOperations, zero_l2_vector_norm)
+{
+    Vector<double> v("v", 3);
+    v(0) = 0;
+    v(1) = 0;
+    v(2) = 0;
+    ConstVector<double> const_v(v);
+    EXPECT_DOUBLE_EQ(l2_norm(ConstVector<double>(const_v)), 0.0);
+}
 
-TEST(VectorOperations, l2_vector_norm)
+TEST(VectorOperations, underflow_l2_vector_norm)
 {
     Vector<double> v("v", 3);
-    v(0) = 1;
-    v(1) = -5;
-    v(2) = 2;
+    v(0) = 1e-300;
+    v(1) = 1e-300;
+    v(2) = 1e-300;
     ConstVector<double> const_v(v);
-    EXPECT_DOUBLE_EQ(l2_norm(ConstVector<double>(const_v)), std::sqrt(30.0));
+    EXPECT_DOUBLE_EQ(l2_norm(ConstVector<double>(const_v)), std::sqrt(3.0) * 1e-300);
+}
+
+TEST(VectorOperations, overflow_l2_vector_norm)
+{
+    Vector<double> v("v", 3);
+    v(0) = 1e+300;
+    v(1) = 1e+300;
+    v(2) = 1e+300;
+    ConstVector<double> const_v(v);
+    EXPECT_DOUBLE_EQ(l2_norm(ConstVector<double>(const_v)), std::sqrt(3.0) * 1e+300);
 }
 
 /* T infinity_norm(ConstVector<T>& x); */