Unitful.jl icon indicating copy to clipboard operation
Unitful.jl copied to clipboard

make assignment work on ustrip'ed StridedArrays

Open bjarthur opened this issue 2 years ago • 10 comments

fixes https://github.com/PainterQubits/Unitful.jl/issues/644

note that i'm only superficially familiar with Unitful.jl, so am not sure if this breaks anything or is consistent with the roadmap. but the tests pass locally.

bjarthur avatar Apr 14 '23 21:04 bjarthur

Codecov Report

Merging #645 (c32177e) into master (4ce89e1) will not change coverage. The diff coverage is 100.00%.

:exclamation: Current head c32177e differs from pull request most recent head 6609979. Consider uploading reports for the commit 6609979 to get more accurate results

:mega: This organization is not using Codecov’s GitHub App Integration. We recommend you install it so Codecov can continue to function properly for your repositories. Learn more

@@           Coverage Diff           @@
##           master     #645   +/-   ##
=======================================
  Coverage   89.07%   89.07%           
=======================================
  Files          16       16           
  Lines        1492     1492           
=======================================
  Hits         1329     1329           
  Misses        163      163           
Impacted Files Coverage Δ
src/utils.jl 95.34% <100.00%> (ø)

:mega: We’re building smart automated test selection to slash your CI/CD build times. Learn more

codecov-commenter avatar Apr 15 '23 06:04 codecov-commenter

This will also apply to ranges, which makes it different from the current behavior:

Before:

julia> r = (1:5)u"m"
(1:5) m

julia> ustrip(r)
1:1:5

After:

julia> r = (1:5)u"m"
(1:5) m

julia> ustrip(r)
5-element reinterpret(Int64, ::StepRange{Quantity{Int64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Quantity{Int64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}):
 1
 2
 3
 4
 5

We could add a method ustrip(r::AbstractRange) = ustrip.(r) to keep the old behavior for ranges, but there are probably more types of AbstractArrays for which reinterpreting isn’t useful, like SArray from the StaticArrays package. Or does it not make a difference in that case (performance-wise)?

Reinterpreting is useful in two ways: to enable mutation of the original array, and to avoid allocations. For ranges and SArrays, both of these do not apply, because they are immutable and stack-allocated (unless you use mutable number types, which you could technically do).

We could implement ustrip(::typeof(reinterpret), ::AbstractArray{<:Quantity}) (see https://github.com/PainterQubits/Unitful.jl/issues/630#issuecomment-1482765400) instead, that way we don’t add to the confusion of when reinterpreting happens and when it doesn’t.

sostock avatar Apr 17 '23 08:04 sostock

how about redefining ustrip(::Array... to dispatch on DenseArray instead of, as this PR currently stands, AbstractArray? doesn't make sense to reinterpret anything other than a dense array anyway, and doing so would preclude StaticArrays and Ranges, the two types you bring up, yet still cover GPUArrays, which is what motivates me.

julia> typeof(SVector(1,2,3)) <: DenseArray
false

julia> typeof(CuVector{Float64}([1,2,3])) <: DenseArray
true

julia> typeof(1:5) <: DenseArray
false

bjarthur avatar Apr 17 '23 15:04 bjarthur

how about redefining ustrip(::Array... to dispatch on DenseArray instead

I like that. Maybe we could even use StridedArray instead of DenseArray, which is more general and still makes sense. BLAS functions usually support StridedArrays.

sostock avatar Apr 18 '23 14:04 sostock

i tried StridedArrays and it works for my use case. have rebased the change into this PR. tests pass locally.

bjarthur avatar Apr 18 '23 23:04 bjarthur

anything else i need to do here? the CI is "awaiting approval".

bjarthur avatar Apr 24 '23 12:04 bjarthur

I just noticed that this is technically a breaking change. At least I consider it as one. A user that calls ustrip(view(…)) in their code may now unintentionally modify the original array, because ustrip used to create a new array before this PR.

This could be avoided if we implement this as ustrip(reinterpret, x). But if we do that, I would want it to error in cases where a wrong result would be returned (#630), which will require some more code.

sostock avatar Apr 24 '23 20:04 sostock

ustrip used to create a new array before this PR because of a deprecated method.

perhaps it's time to bump the version to 2.0? what other new features are being held up because it would be a breaking change that we could implement if we did so?

certainly i'd like to see that deprecated method removed in v2.0. that should've been done in v0.5.

bjarthur avatar Apr 24 '23 21:04 bjarthur

I would like to release v2.0, but there are more breaking changes that I think should be included (and some of them require some further discussion), so it will take some time.

I still think we should use ustrip(reinterpret, x) for this functionality. It seems error-prone to me that ustrip(x) would return an object that shares its memory with x. Maybe I am biased because I did not know that this functionality existed until you opened #627, but I think many users don’t know of this feature, so it has the potential to lead to many errors. I think avoiding errors is more important than convenience.

Functions like view or reshape exist precisely for the purpose of reusing memory and allowing mutation of the original array. But ustrip is mainly a function that is applied to single values, not arrays. Therefore, I think it is not obvious that ustrip(x) might reuse memory.

Is there a reason why you think we shouldn’t add ustrip(reinterpret, x)? If you just think it is too long, maybe we could call it ustrip!(x) instead?

sostock avatar Apr 25 '23 09:04 sostock

should probably add the v2.0 label to this PR or the issue it fixes. i see there are two others with that label.

i don't find it confusing that ustrip on an array reinterprets. that is clearly documented and has existed for a long time. moreover, it is consistent with how arrays are passed in by pointer to functions. and it is probably what one would want to do normally anyway.

what is confusing is that ustriping a unit-less array or a strided one with units allocates.

that we disagree is an indication that more opinions should be solicited. who could we ping that is a heavy user of Unitful.jl?

bjarthur avatar Apr 25 '23 18:04 bjarthur