thchr / Crystalline.jl

Tools for crystallographic symmetry analysis
Other
54 stars 8 forks source link

Compatiblity with Unitful ? #56

Open Joh11 opened 2 months ago

Joh11 commented 2 months ago

Hi,

As we can see with this minimal example, reciprocalbasis does not work well with Unitful quantities:

using Unitful
using Bravais: reciprocalbasis

Rs = ([1.0, 0.0, 0.0], [-0.5, sqrt(3)/2, 0.0],   [0, 0, 1.25])
Us = map(Rs) do v v * 1u"Å" end

@show reciprocalbasis(Rs) # works
@show reciprocalbasis(Us) # StackOverflowError

It yields a stack overflow error. For me, this is likely due to the fact that the input and output dimensions are different ($L^3$ and $L^{-3}$ respectively), so type inference has some troubles.

Is there a possible workaround ? (apart from the trivial one of stripping the units before and putting them back after)

thchr commented 2 months ago

Hmm, the implementation hits a recursive loop because I didn't account for the possibility of the eltypes of the input not being subtypes of Real.

The problematic bits are here: https://github.com/thchr/Crystalline.jl/blob/820b68775ca7076d7cf37a35ad290eadb69d5c3b/Bravais/src/transform.jl#L275 https://github.com/thchr/Crystalline.jl/blob/820b68775ca7076d7cf37a35ad290eadb69d5c3b/Bravais/src/transform.jl#L251-L252

The question, I suppose, is whether DirectBasis and ReciprocalBasis optionally should take a eltype which isn't Float64? Feedback and usecases would be welcome!

thchr commented 2 months ago

So, just to be more specific, the issue here is that

T = eltype(eltype(Us)) # Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}
T <: Real # false

The current implementation assumes all vectors have elements that are subtypes of Real (and Quantity is not).

Related: https://github.com/PainterQubits/Unitful.jl/issues/680