Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

faster and simpler generic_norm2 #43256

Closed
wants to merge 11 commits into from
Closed

faster and simpler generic_norm2 #43256

Changes from 1 commit
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
3e696ed
faster and simpler generic_norm2
oscardssmith Nov 29, 2021
604adc3
accumulate norm2 in at least Float64 precision and convert generic_normp
oscardssmith Dec 2, 2021
44509da
Update stdlib/LinearAlgebra/src/generic.jl
oscardssmith Dec 3, 2021
aedb519
fix typo
oscardssmith Dec 3, 2021
b8c8b85
fix typo
oscardssmith Dec 3, 2021
689032c
Update stdlib/LinearAlgebra/src/generic.jl
oscardssmith Dec 3, 2021
812bd57
fix norm of subnormals
oscardssmith Dec 3, 2021
4306643
test was backwards
oscardssmith Dec 3, 2021
3de392a
fix type stability
oscardssmith Dec 3, 2021
d9de30b
be more lazy
oscardssmith Dec 28, 2021
543f563
Update stdlib/LinearAlgebra/src/generic.jl
oscardssmith Dec 28, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
21 changes: 3 additions & 18 deletions stdlib/LinearAlgebra/src/generic.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -462,27 +462,12 @@ norm_sqr(x::Union{T,Complex{T},Rational{T}}) where {T<:Integer} = abs2(float(x))

function generic_norm2(x)
maxabs = normInf(x)
(maxabs == 0 || isinf(maxabs)) && return maxabs
(v, s) = iterate(x)::Tuple
(ismissing(maxabs) || maxabs == 0 || isinf(maxabs)) && return maxabs
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm sad if we need to add a bunch of special-cased missing code to LinearAlgebra.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I know. We can leave out that test. It's just that the transition to mapreduce facilitated norm handle missing for certain ps, without pulling *missing to the surface.

Copy link
Member

@martinholters martinholters Dec 1, 2021 •
edited by dkarrasch
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Would

Suggested change
(ismissing(maxabs) || maxabs == 0 || isinf(maxabs)) && return maxabs
(isinf(maxabs) !== false || maxabs == 0) && return maxabs

be preferable? (Which I think should be equivalent, but please do double-check.)

EDIT by @dkarrasch : one needs to avoid performing maxabs == 0 because that throws a TypeError: non-boolean (Missing) used in boolean context. Otherwise that would work, because of short-circuiting.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So that you get a notification: I modified @martinholters's suggestion, which I tested for maxabs = missing and it works.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh, and that same trick should make it possible to apply the same steps to generic_normp!

T = typeof(maxabs)
if isfinite(length(x)*maxabs*maxabs) && maxabs*maxabs != 0 # Scaling not necessary
sum::promote_type(Float64, T) = norm_sqr(v)
while true
y = iterate(x, s)
y === nothing && break
(v, s) = y
sum += norm_sqr(v)
end
return convert(T, sqrt(sum))
return convert(T, sqrt(mapreduce(norm_sqr, +, x)))
else
sum = abs2(norm(v)/maxabs)
while true
y = iterate(x, s)
y === nothing && break
(v, s) = y
sum += (norm(v)/maxabs)^2
end
return convert(T, maxabs*sqrt(sum))
return convert(T, maxabs*sqrt(mapreduce(v -> abs2(norm(v)/maxabs), +, x)))
end
end

Expand Down