Handle case when norm==Inf in normalize (#29691)

Fixes #29681

stdlib/LinearAlgebra/src/generic.jl

@@ -1361,17 +1361,18 @@ end
     # The largest positive floating point number whose inverse is less than infinity
     δ = inv(prevfloat(typemax(nrm)))
 
-    if nrm ≥ δ # Safe to multiply with inverse
+    if nrm == Inf  # Just divide since performance isn't important in this corner case
+        v ./= Inf
+    elseif nrm ≥ δ # Safe to multiply with inverse
         invnrm = inv(nrm)
         rmul!(v, invnrm)
-
-    else # scale elements to avoid overflow
+    else           # Scale elements to avoid overflow
         εδ = eps(one(nrm))/δ
         rmul!(v, εδ)
         rmul!(v, inv(nrm*εδ))
     end
 
-    v
+    return v
 end
 
 """

stdlib/LinearAlgebra/test/generic.jl

@@ -256,6 +256,13 @@ end
     @test norm(normalize!(v) - w, Inf) < eps()
 end
 
+@testset "normalize with Infs. Issue 29681." begin
+    @test all(isequal.(normalize([1, -1, Inf]),
+                       [0.0, -0.0, NaN]))
+    @test all(isequal.(normalize([complex(1), complex(0, -1), complex(Inf, -Inf)]),
+                       [0.0 + 0.0im, 0.0 - 0.0im, NaN + NaN*im]))
+end
+
 @testset "Issue 14657" begin
     @test det([true false; false true]) == det(Matrix(1I, 2, 2))
 end