Comments (8)
giordano
commented on May 27, 2024
Uhm, if I define
Base.one(::Type{Measurement{T}}) where T = one(T)
Base.one(::Type{Complex{Measurement{T}}}) where T = one(Complex{T})then some operations are no longer inferred correctly:
julia> @code_warntype(complex(measurement(3)) ^ 2.5)
MethodInstance for ^(::Complex{Measurement{Float64}}, ::Float64)
from ^(z::Complex{T}, p::S) where {T<:Real, S<:Real} @ Base complex.jl:873
Static Parameters
T = Measurement{Float64}
S = Float64
Arguments
#self#::Core.Const(^)
z::Complex{Measurement{Float64}}
p::Float64
Locals
P::Type{Measurement{Float64}}
Body::Union{Complex{Measurement{Float64}}, ComplexF64}
1 ─ (P = Base.promote_type($(Expr(:static_parameter, 1)), $(Expr(:static_parameter, 2))))
│ %2 = Core.apply_type(Base.Complex, P::Core.Const(Measurement{Float64}))::Core.Const(Complex{Measurement{Float64}})
│ %3 = (%2)(z)::Complex{Measurement{Float64}}
│ %4 = (P::Core.Const(Measurement{Float64}))(p)::Core.PartialStruct(Measurement{Float64}, Any[Float64, Core.Const(0.0), Core.Const(0x0000000000000000), Core.Const(Measurements.Derivatives{Float64}())])
│ %5 = (%3 ^ %4)::Union{Complex{Measurement{Float64}}, ComplexF64}
└── return %5
julia> A = [(14 ± 0.1) (23 ± 0.2); (-12 ± 0.3) (29 ± 0.4)];
julia> using LinearAlgebra
julia> @code_warntype det(A)
MethodInstance for LinearAlgebra.det(::Matrix{Measurement{Float64}})
from det(A::AbstractMatrix{T}) where T @ LinearAlgebra ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/LinearAlgebra/src/generic.jl:1617
Static Parameters
T = Measurement{Float64}
Arguments
#self#::Core.Const(LinearAlgebra.det)
A::Matrix{Measurement{Float64}}
Locals
S::Type{Measurement{Float64}}
Body::Union{Float64, Measurement{Float64}}
1 ─ Core.NewvarNode(:(S))
│ %2 = LinearAlgebra.istriu(A)::Bool
└── goto #3 if not %2
2 ─ goto #4
3 ─ %5 = LinearAlgebra.istril(A)::Bool
└── goto #5 if not %5
4 ┄ %7 = $(Expr(:static_parameter, 1))::Core.Const(Measurement{Float64})
│ %8 = LinearAlgebra.one($(Expr(:static_parameter, 1)))::Core.Const(1.0)
│ %9 = LinearAlgebra.zero($(Expr(:static_parameter, 1)))::Core.Const(0.0 ± 0.0)
│ %10 = (%8 * %9)::Core.Const(0.0 ± 0.0)
│ %11 = LinearAlgebra.zero($(Expr(:static_parameter, 1)))::Core.Const(0.0 ± 0.0)
│ %12 = (%10 + %11)::Core.PartialStruct(Measurement{Float64}, Any[Core.Const(0.0), Float64, Core.Const(0x0000000000000000), Measurements.Derivatives{Float64}])
│ %13 = LinearAlgebra.one($(Expr(:static_parameter, 1)))::Core.Const(1.0)
│ %14 = (%12 / %13)::Core.PartialStruct(Measurement{Float64}, Any[Core.Const(0.0), Float64, Core.Const(0x0000000000000000), Measurements.Derivatives{Float64}])
│ %15 = LinearAlgebra.typeof(%14)::Core.Const(Measurement{Float64})
│ (S = LinearAlgebra.promote_type(%7, %15))
│ %17 = S::Core.Const(Measurement{Float64})
│ %18 = LinearAlgebra.UpperTriangular(A)::UpperTriangular{Measurement{Float64}, Matrix{Measurement{Float64}}}
│ %19 = LinearAlgebra.det(%18)::Core.PartialStruct(Measurement{Float64}, Any[Float64, Float64, Core.Const(0x0000000000000000), Measurements.Derivatives{Float64}])
│ %20 = LinearAlgebra.convert(%17, %19)::Core.PartialStruct(Measurement{Float64}, Any[Float64, Float64, Core.Const(0x0000000000000000), Measurements.Derivatives{Float64}])
└── return %20
5 ─ %22 = (:check,)::Core.Const((:check,))
│ %23 = Core.apply_type(Core.NamedTuple, %22)::Core.Const(NamedTuple{(:check,)})
│ %24 = Core.tuple(false)::Core.Const((false,))
│ %25 = (%23)(%24)::Core.Const((check = false,))
│ %26 = Core.kwcall(%25, LinearAlgebra.lu, A)::LU{Measurement{Float64}, Matrix{Measurement{Float64}}, Vector{Int64}}
│ %27 = LinearAlgebra.det(%26)::Union{Float64, Measurement{Float64}}
└── return %27
which breaks tests at
Measurements.jl/test/runtests.jl
Lines 288 to 289 in daa1ad2
Measurements.jl/test/runtests.jl
Line 732 in daa1ad2
from measurements.jl.
giordano
commented on May 27, 2024
And the accuracy of a bunch of QuadGK tests gets worse:
QuadGK: Test Failed at /home/mose/.julia/dev/Measurements/test/runtests.jl:852
Expression: ≈((QuadGK.quadgk(sin, -y, y))[1], cos(-y) - cos(y), atol = eps(Float64))
Evaluated: -2.23e-16 ± 2.7e-17 ≈ 0.0 ± 0.0 (atol=2.220446049250313e-16)
Stacktrace:
[1] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:477 [inlined]
[2] macro expansion
@ ~/.julia/dev/Measurements/test/runtests.jl:852 [inlined]
[3] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:1496 [inlined]
[4] top-level scope
@ ~/.julia/dev/Measurements/test/runtests.jl:851
QuadGK: Test Failed at /home/mose/.julia/dev/Measurements/test/runtests.jl:855
Expression: ≈((QuadGK.quadgk((t->(cos(x - t);)), 0, 2pi))[1], measurement(0), atol = 7.0e-16)
Evaluated: -1.063e-13 ± 2.9e-15 ≈ 0.0 ± 0.0 (atol=7.0e-16)
Stacktrace:
[1] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:477 [inlined]
[2] macro expansion
@ ~/.julia/dev/Measurements/test/runtests.jl:855 [inlined]
[3] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:1496 [inlined]
[4] top-level scope
@ ~/.julia/dev/Measurements/test/runtests.jl:851
QuadGK: Test Failed at /home/mose/.julia/dev/Measurements/test/runtests.jl:869
Expression: ≈((QuadGK.quadgk(sin, -y, y))[1], #= /home/mose/.julia/dev/Measurements/test/runtests.jl:870 =# @uncertain(((x->begin
(QuadGK.quadgk(sin, -x, x))[1]
end))(y)), atol = 1.0e-10)
Evaluated: -2.23e-16 ± 2.7e-17 ≈ -2.2e-16 ± 4.6e-10 (atol=1.0e-10)
Stacktrace:
[1] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:477 [inlined]
[2] macro expansion
@ ~/.julia/dev/Measurements/test/runtests.jl:869 [inlined]
[3] macro expansion
@ ~/.julia/juliaup/julia-1.9.0-beta4+0.x64.linux.gnu/share/julia/stdlib/v1.9/Test/src/Test.jl:1496 [inlined]
[4] top-level scope
@ ~/.julia/dev/Measurements/test/runtests.jl:851
Test Summary: | Pass Fail Total Time
QuadGK | 16 3 19 14.0s
although the prospect of deleting the src/quadgk.jl file is quite appealing.
from measurements.jl.
giordano
commented on May 27, 2024
I presume the problem with the complex power is the use of one at https://github.com/JuliaLang/julia/blob/afeda9f8cf0776461eaf565333e2159122d70bc4/base/complex.jl#L789?
Edit: no, or not only that: replacing one with oneunit on that line doesn't fix the inference failure.
from measurements.jl.
giordano
commented on May 27, 2024
And the accuracy of a bunch of QuadGK tests gets worse:
Actually, that's independent from this change, accuracy got worse in v2.8.1 of QuadGK, it's restored if I downgraded it to v2.8.0.
from measurements.jl.
giordano
commented on May 27, 2024
In particular, with 2.8.0:
julia> QuadGK.quadgk(t -> cos(3 - t), 0, 2pi)
(6.106226635438361e-16, 0.0)with 2.8.1:
julia> QuadGK.quadgk(t -> cos(3 - t), 0, 2pi)
(-1.0625509660110193e-13, 1.1749641421160047e-23)This has nothing to do with Measurements.jl
from measurements.jl.
giordano
commented on May 27, 2024
@stevengj Additionally, v2.6.0 of QuadGK broke another test here:
Measurements.jl/test/runtests.jl
Lines 925 to 930 in daa1ad2
julia> using Unitful, Measurements, QuadGK
julia> f(x) = x^2
f (generic function with 1 method)
julia> F(x) = x^3 / 3
F (generic function with 1 method)
julia> a = (5 ± 0.1)u"m"
5.0 ± 0.1 m
julia> b = (10 ± 1)u"m"
10.0 ± 1.0 m
julia> QuadGK.quadgk(f, a, b)
ERROR: MethodError: no method matching Float64(::Measurement{Float64})
Closest candidates are:
(::Type{T})(::Real, ::RoundingMode) where T<:AbstractFloat at rounding.jl:200
(::Type{T})(::T) where T<:Number at boot.jl:772
(::Type{T})(::AbstractChar) where T<:Union{AbstractChar, Number} at char.jl:50
...
Stacktrace:
[1] convert(#unused#::Type{Float64}, x::Measurement{Float64})
@ Base ./number.jl:7
[2] setindex!(A::Vector{Float64}, x::Measurement{Float64}, i1::Int64)
@ Base ./array.jl:966
[3] eignewt(b::Vector{Measurement{Float64}}, m::Int64, n::Int64)
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/gausskronrod.jl:43
[4] kronrod(#unused#::Type{Measurement{Float64}}, n::Int64)
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/gausskronrod.jl:197
[5] macro expansion
@ ~/.julia/packages/QuadGK/gxzkm/src/gausskronrod.jl:259 [inlined]
[6] _cachedrule(#unused#::Type{Measurement{Float64}}, n::Int64)
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/gausskronrod.jl:259
[7] cachedrule
@ ~/.julia/packages/QuadGK/gxzkm/src/gausskronrod.jl:264 [inlined]
[8] do_quadgk(f::typeof(f), s::Tuple{Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}, n::Int64, atol::Nothing, rtol::Nothing, maxevals::Int64, nrm::typeof(LinearAlgebra.norm), segbuf::Nothing)
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/adapt.jl:7
[9] #46
@ ~/.julia/packages/QuadGK/gxzkm/src/adapt.jl:219 [inlined]
[10] handle_infinities(workfunc::QuadGK.var"#46#47"{Nothing, Nothing, Int64, Int64, typeof(LinearAlgebra.norm), Nothing}, f::typeof(f), s::Tuple{Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}})
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/adapt.jl:118
[11] quadgk(::Function, ::Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, ::Vararg{Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}; atol::Nothing, rtol::Nothing, maxevals::Int64, order::Int64, norm::Function, segbuf::Nothing)
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/adapt.jl:218
[12] quadgk(::Function, ::Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, ::Quantity{Measurement{Float64}, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}})
@ QuadGK ~/.julia/packages/QuadGK/BYxcx/src/adapt.jl:216
[13] top-level scope
@ REPL[36]:1Although this would be fixed in #135 by defining the one method you suggested.
from measurements.jl.
stevengj
commented on May 27, 2024
QuadGK.quadgk(t -> cos(3 - t), 0, 2pi)
The basic problem with this test is that you are asking it to compute an integral that is zero to a relative tolerance of ≈ 1e-8 (the default), which is basically impossible. So, the convergence criterion in practice depends on arbitrary details of roundoff errors, so e.g. slightly changing the order of the integrand evaluation (as we did in JuliaMath/QuadGK.jl@962c801 for 2.8.1) causes it to stop at a very different point.
The solution, as explained in the QuadGK manual, is that whenever you have an integral that might be zero, you should pass in a nonzero atol. Then it works fine in both 2.8.1, though it still can't converge to better than machine precision:
julia> quadgk_count(t -> cos(3 - t), 0, 2pi, atol=1e-15)
(2.220446049250313e-16, 5.551115123125783e-16, 105)
and 2.8.0:
julia> quadgk_count(t -> cos(3 - t), 0, 2pi, atol=1e-15)
(3.3306690738754696e-16, 2.220446049250313e-16, 105)
though the answers are very slightly different due to changes in roundoff.
from measurements.jl.
stevengj
commented on May 27, 2024
For powers, oneunit is not correct since if you have a unitful number x then x^0 should be unitless, i.e. the x^n function is inherently type-unstable for dimensionful numbers.
The alternative would be to define some other method to get the "scalar" type from a numeric type. Maybe call it scalarone(x) or something like that, then define a package ScalarOne.jl and try to get all of the packages (Unitful, ForwardDiff, etcetera) to adopt it. It could call scalarone(x::Number) = one(x) as a fallback, so it would be no worse than what we have now, I guess.
from measurements.jl.
Related Issues (20)
- Adding measurement components back to a measurement after iteratively solving for a value HOT 7
- [FR] Plot recipe: ribbon plots option beside error bars HOT 11
- tryparse for Measurement type
- Can't use unique with measurements HOT 1
- Measurements with missing errors HOT 4
- Measurements.value(x::Missing) = missing HOT 1
- Integration with Zygote
- `weightedmean()` returns `NaN ± 0`? HOT 1
- Use auto-differentiation engine
- Bad integration with Plots' boxplot HOT 2
- Move to pkgextensions for Julia v1.9+
- `Symbolics.jl` support? HOT 7
- Plotting mixture of measurements and missing data HOT 1
- Trying to use Measurements to differentiate respect to a unitful quantity. HOT 12
- Is there an autodiff package that is compatible? HOT 3
- Broken `MeasurementsJunoExt.jl` HOT 1
- Julia 1.6 incompatibility from stdlib compat bounds HOT 2
- Wrong values when multiplication with scalar? HOT 3
- Distinguishing between statistical and systematic errors in Measurements.jl HOT 1
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from measurements.jl.