diff --git a/README.md b/README.md index f2d7685a..0ad9e142 100644 --- a/README.md +++ b/README.md @@ -326,8 +326,8 @@ julia> ustrip(ua"degC", 295.15u"K") For working with an array of quantities that have the same dimensions, you can use a `QuantityArray`: -```julia -julia> ar = QuantityArray(rand(3), u"m/s") +```julia-repl +julia> ar = rand(3)u"m/s" 3-element QuantityArray(::Vector{Float64}, ::Quantity{Float64, Dimensions{FixedRational{Int32, 25200}}}): 0.2729202669351497 m s⁻¹ 0.992546340360901 m s⁻¹ @@ -352,7 +352,7 @@ This also has a custom broadcasting interface which allows the compiler to avoid redundant dimension calculations, relative to if you had simply used an array of quantities: -```julia +```julia-repl julia> f(v) = v^2 * 1.5; julia> @btime $f.(xa) setup=(xa = randn(100000) .* u"km/s"); @@ -364,6 +364,16 @@ julia> @btime $f.(qa) setup=(xa = randn(100000) .* u"km/s"; qa = QuantityArray(x So we can see the `QuantityArray` version saves on both time and memory. +By default, DynamicQuantities will create a `QuantityArray` from an `AbstractArray`, similarly to how a `Quantity` is created from a scalar in the [Usage](@ref) examples: + +```julia-repl +julia> x = (1:3)us"km/h" +3-element QuantityArray(::StepRangeLen{Float64, ...}, ::Quantity{Float64, SymbolicDimensions{...}}): + 1.0 km h⁻¹ + 2.0 km h⁻¹ + 3.0 km h⁻¹ +``` + ### Unitful DynamicQuantities allows you to convert back and forth from Unitful.jl: diff --git a/docs/src/examples.md b/docs/src/examples.md index e5dc7d1f..19b5a72c 100644 --- a/docs/src/examples.md +++ b/docs/src/examples.md @@ -238,6 +238,7 @@ the same dimension) by passing an array and a single quantity: ```julia x = QuantityArray(randn(32), u"km/s") +# or x = randn(32)u"km/s" ``` or, by passing an array of individual quantities: @@ -273,6 +274,8 @@ f_square(v) = v^2 * 1.5 - v^2 println("Applying function to y_q: ", sum(f_square.(y_q))) ``` +See [Home > Arrays](index.md#Arrays) for more. + ### Fill We can also make `QuantityArray` using `fill`: diff --git a/ext/DynamicQuantitiesLinearAlgebraExt.jl b/ext/DynamicQuantitiesLinearAlgebraExt.jl index 1b2b5699..0cfa4bec 100644 --- a/ext/DynamicQuantitiesLinearAlgebraExt.jl +++ b/ext/DynamicQuantitiesLinearAlgebraExt.jl @@ -2,7 +2,7 @@ module DynamicQuantitiesLinearAlgebraExt using LinearAlgebra: LinearAlgebra as LA using DynamicQuantities -using DynamicQuantities: DynamicQuantities as DQ, quantity_type, new_quantity, DimensionError +using DynamicQuantities: DynamicQuantities as DQ, quantity_type, new_quantity, DimensionError, ABSTRACT_QUANTITY_TYPES using TestItems: @testitem DQ.is_ext_loaded(::Val{:LinearAlgebra}) = true @@ -35,13 +35,37 @@ for op in (:(Base.:*), :(Base.:/), :(Base.:\)), @eval $op(l::$L, r::$R) = DQ.array_op($op, l, r) end +for ARRAY_TYPE in ( + LA.Bidiagonal, + LA.Diagonal, + LA.Hermitian, + LA.LowerTriangular, + LA.LowerTriangular{<:Any, <:Union{LA.Adjoint{<:Any, <:StridedMatrix{T}}, LA.Transpose{<:Any, <:StridedMatrix{T}}, StridedArray{T, 2}} where T}, + LA.Symmetric, + LA.SymTridiagonal, + LA.Tridiagonal, + LA.UnitLowerTriangular, + LA.UnitUpperTriangular, + LA.UpperTriangular, + LA.UpperTriangular{<:Any, <:Union{LA.Adjoint{<:Any, <:StridedMatrix{T}}, LA.Transpose{<:Any, <:StridedMatrix{T}}, StridedArray{T, 2}} where T}, + LA.UpperHessenberg, + ), + (type, _, _) in ABSTRACT_QUANTITY_TYPES + + @eval begin + Base.:*(A::$ARRAY_TYPE, q::$type) = QuantityArray(A, q) + Base.:*(q::$type, A::$ARRAY_TYPE) = A * q + Base.:/(A::$ARRAY_TYPE, q::$type) = A * inv(q) + end +end + function Base.:*( l::LA.Transpose{Q,<:AbstractVector}, r::DQ.QuantityArray{T2,1,D,Q,<:AbstractVector{T2}} ) where { T2,D<:DQ.AbstractDimensions,Q<:DQ.AbstractRealQuantity{T2,D} } - return array_op(Base.:*, l, r) + return DQ.array_op(Base.:*, l, r) end @@ -52,6 +76,30 @@ end end # TODO: functions on SVD type that are working: `size`, `adjoint`, partially working: `inv`, not working: `svdvals`, `ldiv!`. +@testitem "QuantityArray construction" begin + using DynamicQuantities, LinearAlgebra + + tname = nameof ∘ typeof + A = [1 0; 0 1] + @testset "$(tname(arr)) and $(tname(q))" for arr in ( + Bidiagonal(A, :U), + Diagonal(A), + Hermitian(A), + LowerTriangular(A), + Symmetric(A), + SymTridiagonal(A), + Tridiagonal(A), + UnitLowerTriangular(A), + UnitUpperTriangular(A), + UpperTriangular(A), + UpperHessenberg(A), + ), + q in (u"m", GenericQuantity(1.5), RealQuantity(1.5)) + @test arr * q == QuantityArray(A, q) + @test q * arr == QuantityArray(A, q) + @test arr / q == arr * inv(q) + end +end @testitem "svd" begin using DynamicQuantities, LinearAlgebra @@ -211,4 +259,14 @@ end # TODO: Tests from missing parts of LinearAlgebra interface +@testitem "transpose" begin + using DynamicQuantities, LinearAlgebra + + v = QuantityArray([1, 2, 3], RealQuantity(u"m")) + square_norm = transpose(v) * v + + @test square_norm == QuantityArray([14], RealQuantity(u"m^2")) + @test length(square_norm) == 1 +end + end diff --git a/src/arrays.jl b/src/arrays.jl index 53a1522d..a9d0ee1f 100644 --- a/src/arrays.jl +++ b/src/arrays.jl @@ -14,6 +14,14 @@ and so can be used in most places where a normal array would be used, including # Constructors +You can create a `QuantityArray` by multiplying an array with a `Quantity`: + +```julia +x = (0:0.1:10)u"km/s" +``` + +Alternatively, the following constructors are available: + - `QuantityArray(v::AbstractArray, d::AbstractDimensions)`: Create a `QuantityArray` with value `v` and dimensions `d`, using `Quantity` if the eltype of `v` is numeric, and `GenericQuantity` otherwise. @@ -71,6 +79,20 @@ for (type, base_type, default_type) in ABSTRACT_QUANTITY_TYPES if type in (AbstractQuantity, AbstractGenericQuantity) @eval QuantityArray(v::AbstractArray{<:$base_type}, d::AbstractDimensions) = QuantityArray(v, d, $default_type) end + + @eval begin + Base.:*(A::AbstractArray{T}, q::$type) where {T<:Number} = QuantityArray(A, q) + Base.:*(q::$type, A::AbstractArray{T}) where {T<:Number} = A * q + Base.:/(A::AbstractArray{T}, q::$type) where {T<:Number} = A * inv(q) + end + + for (type2, _, _) in ABSTRACT_QUANTITY_TYPES + @eval begin + Base.:*(A::AbstractArray{<:$type2}, q::$type) = A .* q + Base.:*(q::$type, A::AbstractArray{<:$type2}) = q .* A + Base.:/(A::AbstractArray{<:$type2}, q::$type) = A ./ q + end + end end QuantityArray(v::QA) where {Q<:UnionAbstractQuantity,QA<:AbstractArray{Q}} = let diff --git a/src/disambiguities.jl b/src/disambiguities.jl index 7896c555..7d822b91 100644 --- a/src/disambiguities.jl +++ b/src/disambiguities.jl @@ -110,3 +110,9 @@ for (type, _, _) in ABSTRACT_QUANTITY_TYPES, numeric_type in (Bool, BigFloat) end end end +# Cover method ambiguities from, e.g., op(::Array, ::Quantity)::QuantityArray` +Base.:*(A::StepRangeLen{<:Real, <:Base.TwicePrecision}, q::AbstractRealQuantity) = QuantityArray(A, q) +Base.:*(q::AbstractRealQuantity, A::StepRangeLen{<:Real, <:Base.TwicePrecision}) = A * q +Base.:/(A::BitArray, q::AbstractRealQuantity) = A * inv(q) +Base.:/(A::BitArray, q::AbstractQuantity) = A * inv(q) +Base.:/(A::StepRangeLen{<:Real, <:Base.TwicePrecision}, q::AbstractRealQuantity) = A * inv(q) diff --git a/test/unittests.jl b/test/unittests.jl index 189dab09..d0d9dc2f 100644 --- a/test/unittests.jl +++ b/test/unittests.jl @@ -288,9 +288,9 @@ end @test ustrip(x + ones(T, 32))[32] == 2 @test typeof(x + ones(T, 32)) <: GenericQuantity{Vector{T}} @test typeof(x - ones(T, 32)) <: GenericQuantity{Vector{T}} - @test typeof(ones(T, 32) * GenericQuantity(T(1), D, length=1)) <: GenericQuantity{Vector{T}} - @test typeof(ones(T, 32) / GenericQuantity(T(1), D, length=1)) <: GenericQuantity{Vector{T}} - @test ones(T, 32) / GenericQuantity(T(1), length=1) == GenericQuantity(ones(T, 32), length=-1) + @test typeof(ones(T, 32) * GenericQuantity(T(1), D, length=1)) <: QuantityArray + @test typeof(ones(T, 32) / GenericQuantity(T(1), D, length=1)) <: QuantityArray + @test ones(T, 32) / GenericQuantity(T(1), length=1) == QuantityArray(ones(T, 32), GenericQuantity(T(1), length=-1)) end @testset "isapprox" begin @@ -351,6 +351,11 @@ end @test norm(GenericQuantity(ustrip.(x), length=1, time=-1), 2) ≈ norm(ustrip.(x), 2) * u"m/s" @test ustrip(x') == ustrip(x)' + + # With BitArray and RealQuantity: + @test BitArray([1 0 1 0]) / u"m" isa QuantityArray + @test BitArray([1 0 1 0]) / us"m" isa QuantityArray + @test BitArray([1 0 1 0]) / RealQuantity(u"m") isa QuantityArray end @testset "Ranges" begin @@ -397,7 +402,7 @@ end @testset "Multiplying ranges with units" begin # Test multiplying ranges with units x = (1:0.25:4)u"inch" - @test x isa StepRangeLen + @test x isa QuantityArray @test first(x) == 1u"inch" @test x[2] == 1.25u"inch" @test last(x) == 4u"inch" @@ -405,7 +410,7 @@ end # Integer range (but real-valued unit) x = (1:4)u"inch" - @test x isa StepRangeLen + @test x isa QuantityArray @test first(x) == 1u"inch" @test x[2] == 2u"inch" @test last(x) == 4u"inch" @@ -414,7 +419,7 @@ end # Test with floating point range x = (1.0:0.5:3.0)u"m" - @test x isa StepRangeLen + @test x isa QuantityArray @test first(x) == 1.0u"m" @test x[2] == 1.5u"m" @test last(x) == 3.0u"m" @@ -427,7 +432,7 @@ end # Test with symbolic units x = (1:0.25:4)us"inch" - @test x isa StepRangeLen{<:Quantity{Float64,<:SymbolicDimensions}} + @test x isa QuantityArray @test first(x) == us"inch" @test x[2] == 1.25us"inch" @test last(x) == 4us"inch" @@ -435,7 +440,7 @@ end # Test that symbolic units preserve their symbolic nature x = (0:0.1:1)us"km/h" - @test x isa AbstractRange + @test x isa QuantityArray @test first(x) == 0us"km/h" @test x[2] == 0.1us"km/h" @test last(x) == 1us"km/h" @@ -443,7 +448,7 @@ end # Similarly, integers should stay integers: x = (1:4)us"inch" - @test_skip x isa StepRangeLen{<:Quantity{Int64,<:SymbolicDimensions}} + @test_skip x isa QuantityArray @test first(x) == us"inch" @test x[2] == 2us"inch" @test last(x) == 4us"inch" @@ -453,6 +458,44 @@ end @test_skip (1.0:4.0) * RealQuantity(u"inch") isa StepRangeLen{<:RealQuantity{Float64,<:SymbolicDimensions}} # TODO: This is not available as TwicePrecision interacts with Real in a way # that demands many other functions to be defined. + @test (1.0:4.0) / RealQuantity(u"inch") isa QuantityArray + @test (1.0:4.0) * RealQuantity(u"inch") isa QuantityArray + @test RealQuantity(u"inch") * (1.0:4.0) isa QuantityArray + + @test (1.0:4.0) / RealQuantity(us"inch") isa QuantityArray + @test (1.0:4.0) * RealQuantity(us"inch") isa QuantityArray + @test RealQuantity(us"inch") * (1.0:4.0) isa QuantityArray + end + + @testset "Array of quantities multiplication behavior" begin + # Test that array of quantities * quantity does NOT create nested QuantityArray + # but instead broadcasts properly + + # Test all combinations of quantity types + for (Q1, Q2) in [(Quantity, Quantity), (Quantity, RealQuantity), + (RealQuantity, Quantity), (RealQuantity, RealQuantity), + (GenericQuantity, GenericQuantity), (GenericQuantity, Quantity)] + + # Create array of Q1 quantities + arr = [Q1(1.0u"m"), Q1(2.0u"m"), Q1(3.0u"m")] + q = Q2(1.0u"s") + + # Test array * quantity + result1 = arr * q + @test result1 isa Vector # Should be Vector, not QuantityArray + @test length(result1) == 3 + @test dimension(result1[1]) == dimension(u"m*s") + + # Test quantity * array + result2 = q * arr + @test result2 isa Vector + @test dimension(result2[1]) == dimension(u"m*s") + + # Test division + result3 = arr / q + @test result3 isa Vector + @test dimension(result3[1]) == dimension(u"m/s") + end end end