Difficulty defining promotion rules

General Usage
1

Hi,

I’m trying to define promotion rules for a type representing multicomplex numbers (Multicomplex number - Wikipedia). These have a series of complex numbers, im1, im2, etc. Internally, these are represented as static arrays of size 2^N:

using StaticArrays
struct Multicomplex{T,N,C} <: Number
    value::SVector{C,T}
end
im1 = Multicomplex{Bool,1,2}(SVector{2}(false, true))
im2 = Multicomplex{Bool,2,4}(SVector{4}(false, false, true, false))

I’m trying to define a promotion rule so that Multicomplex numbers (e.g. im1) get promoted to larger sizes when needed in calculations:

promote_rule(::Type{<:Multicomplex{T,N1,C1}},
             ::Type{<:Multicomplex{S,N2,C2}}) where
             {T<:Real,S<:Real,N1,C1,N2,C2} =
     Multicomplex{promote_type(T,S),maximum(N1,N2),maximum(C1,C2)}
# constructors are also defined for conversion, e.g.
#Multicomplex{T,2,4}(m::Multicomplex{S,1,2}) where {T,S} =
#    Multicomplex{2}(SVector{4,T}(m.value[1], m.value[2], zero(T), zero(T)))

However, this rule doesn’t seem to work, as type parameters are just dropped when they don’t match in promote_type:

julia> promote_type(Multicomplex{Float64,1,2}, Multicomplex{Float64,1,2})
Multicomplex{Float64, 1, 2}

julia> promote_type(Multicomplex{Float64,1,2}, Multicomplex{Float64,1,4})
Multicomplex{Float64, 1}

julia> promote_type(Multicomplex{Float64,1,2}, Multicomplex{Float64,2,4})
Multicomplex{Float64}

julia> promote_type(Multicomplex{Float64,1,2}, Multicomplex{Bool,2,4})
Multicomplex

Any advice would be appreciated. Thank you!

2

There are two problems there.

First, you need to import Base.promote_rule (your definition is being ignored by promote_type). Second, maximum is for arrays, you want max:


julia> function Base.promote_rule(::Type{<:Multicomplex{T,N1,C1}},
                    ::Type{<:Multicomplex{S,N2,C2}}) where
                    {T<:Real,S<:Real,N1,C1,N2,C2} 
            Multicomplex{promote_type(T,S),max(N1,N2),max(C1,C2)}
       end

julia> promote_type(Multicomplex{Float64,1,2}, Multicomplex{Float64,1,4})
Multicomplex{Float64, 1, 4}
3

Two stupid mistakes :grinning: I always get max and maximum confused, and had imported most other functions from Base but forgot this one! Thank you for the quick response!

4

What does the second type parameter, N, encode? It seems redundant if it’s just C = 2^N. Can’t you use only the C parameter for simplicity?

Another thing, did you consider using Tuple instead of SVector? Or do you need StaticArray functionality?

5

Yes, it is encoding C=2^N. My thinking is this: I’d originally tried to avoid multiple parameters by having type {T,N} only, but then I couldn’t correctly set the size of the SVector. I suppose it’s possible to have the type as {T,C}, but most users will be thinking in terms of N and not C, so I’d rather expose this. By having the type as {T,N,C} (and constructors specified in terms of T and N), users can ignore the last parameter C and just work with Multicomplex{T,N} objects.

If you’ve a better suggestion though it would be welcome - I’m a relative newbie.

My thinking for SVector was that I’ll be doing arithmetic on the contents, e.g.

Base.:(+)(a::Multicomplex{T,N}, b::Multicomplex{S,N}) where {T,S,N} =
    Multicomplex{N}(a.value + b.value)
Base.:(-)(a::Multicomplex{T,N}, b::Multicomplex{S,N}) where {T,S,N} =
    Multicomplex{N}(a.value - b.value)

This wouldn’t be so straightforward with Tuple?

6

I see. It’s a bit hard to say, but it seemed like values would be constructed as

im2 = Multicomplex{Bool,2,4}(SVector{4}(false, false, true, false))

etc, so then the user would anyway have to consider the number of values to input, and think ‘in terms of C’. I guess it depends on how this will be used in the end.

Tuples can also do arithmetic with the help of broadcasting

julia> (1,2,3) .+ (4,5,6)
(5, 7, 9)

which should be as performant as staticarrays (StaticArrays.jl is built on tuples.)

One problem, at least with the code you shared, is that nothing is stopping me from writing

im1 = Multicomplex{Bool,11,2}(SVector{2}(false, true))
im2 = Multicomplex{Bool,47,2}(SVector{2}(false, true))
# etc
7

Multicomplex numbers are often built-up (and mathematically defined) recursively from pairs of lower-order numbers:

m_(n) = a_(n-1) + i_n b_(n-1)

so this route to construction avoids thinking of the number of components (directly, at least). It’s a grey area though, and I agree that there’s definitely times when it’s hard to avoid thinking about the number of components. I do feel it’s generally going to be preferable to expose N to the user rather than C though.

I’ve got some checks in the constructors that I didn’t show to avoid invalid pairs of N and C - and you’d still need similar checks to ensure C was a power of two.

That’s helpful to know about Tuples and arithmetic - it would be good to reduce dependence on other packages. One other area I use StaticArrays is when gradually constructing vectors via an MVector, e.g. to pad a smaller vector with zeros:

v = zeros(MVector{C,T})
for i=1:C0
    v[i] = m.value[i]
end
Multicomplex{N}(SVector(v))

I don’t know if there’d be a cleaner way of doing this, particularly with Tuple?

Thanks again - I’m appreciating the discussion!

8

I can see why that’s more aesthetically pleasing, yes.

As for this:

you could try something like this (I mean, I’m not certain what the interface should be, but this is equivalent, I think):

function padto(t::NTuple{N, T}, ::Val{P}) where {N, P, T}
    m = P - N
    return (t..., ntuple(_->zero(T), m)...)
end

I was actually not sure if the compiler would be able to optimize away the P - M, but I think it did. So this is if you want to pad to a certain length. If you want to pad with a certain number of zeros, it’s easier:

function padwith(t::NTuple{N, T}, ::Val{P}) where {N, P, T}
    return (t..., ntuple(_->zero(T), P)...)
end

(You call them as padto(t, Val(7)) to pad to length 7, etc.)

9

Great - that looks really useful. I’ll definitely have a think about switching across. Thanks for all the help!