Set array dimensions in composite type

General Usage
1

How can you define arrays of increasing dimension in a composite type when they have a simple form such as the following?

struct compound{N} # N is the number of spatial dimensions of the field
    scalar_field :: Array{Float64,N}
    vector_field :: Array{Float64,N+1}
    tensor_field :: Array{Float64,N+2}
end

But this throws the error

MethodError: no method matching +(::TypeVar, ::Int64)

I guess I could make an Array of Arrays, but I would rather keep them flat like this.

2

I think you have to settle for three type parameters and an assertion:

struct Compound{N, M, P}
    scalar_field::Array{Float64, N}
    vector_field::Array{Float64, M}
    tensor_field::Array{Float64, P}

    function Compound(s, v, t)
        # alternatively throw a more appropriate error type like
        # a DimensionMismatch or MethodError
        @assert ndims(s) == ndims(v)-1 == ndims(t)-2 "error message about dimension" 
        N = ndims(s)
        new{N, N+1, N+2}(s, v, t)
    end
end
2 Likes
3

I don’t know of any better way than parameterizing Compound with three independent parameters (e.g. N, N1, N2) and enforce the relation between them in an inner constructor.

Something along the lines of:

julia> struct Compound{N, N1, N2}
           scalar_field :: Array{Float64,N}  # N is the number of spatial dimensions of the field
           vector_field :: Array{Float64,N1} # N1 should be N+1
           tensor_field :: Array{Float64,N2} # N2 should be N+2
       
           function Compound(N)
               new{N, N+1, N+2}()
           end
       end

julia> Compound(5)
Compound{5,6,7}(#undef, #undef, #undef)
1 Like
4

Drats! That’s what I’m doing now and it is soooo ugly.

5 
struct Compound{N}
    scalar_field :: Array{Float64,N}
    vector_field :: Array{SVector{1,Float},N}
    tensor_field :: Array{SVector{2,Float},N} 
end

Maybe ?

6

Yeah, the is the array of arrays that I was considering above. Is there a performance hit for nesting arrays like this?

7

You could try with meta programming


abstract type AbstractCompound{N} end

for N = 1:4
    name = Symbol("Compound$N")
    @eval begin
        struct $name <: AbstractCompound{$(N)} # N is the number of spatial dimensions of the field
            scalar_field :: Array{Float64,$(N)}
            vector_field :: Array{Float64,$(N+1)}
            tensor_field :: Array{Float64,$(N+2)}
        end
    end
end

julia> c = Compound1(zeros(1),zeros(1,1),zeros(1,1,1))
Compound1([0.0], [0.0], [0.0])

julia> c.tensor_field
1Ă—1Ă—1 Array{Float64, 3}:
[:, :, 1] =
 0.0
1 Like
8

Oops ! Too fast reading (although not fast enough to provide the first solution :wink: )
I think that you can reshape your arrays as flat arrays.
If you want to minimize the performance risks maybe you may adopt a SoA strategy (but the client code is likely to be “uglier”)

I am thinking about the best strategy and I confess that at this point nothing’s really convince me…

1 Like
9

Neat! This would be much nicer when writing thing like foo(a::AbstractCompound{N}) where N.

Would it be possible to overload the internal constructor? I often have code like b = Compound(a.scalar_field) and now it looks like I would need to check if I needed Compound2, Compound3,…

10

Thanks. That’s my current approach.

11

You can solve this problem the same way, with meta programming. The only thing to keep in mind is to interpolate the N variable into the evaluated expression with $(N).

2 Likes
12

This seems to work well!

abstract type AbstractCompound{N} end
for N = 1:4
    name = Symbol("Compound$N")
    @eval begin
        struct $name <: AbstractCompound{$(N)} # N is the number of spatial dimensions of the field
            scalar_field :: Array{Float64,$(N)}
            vector_field :: Array{Float64,$(N+1)}
            tensor_field :: Array{Float64,$(N+2)}
        end
        Compound(s::AbstractArray{Float64,$(N)},a...) = $name(s,a...)
    end
end
Compound(N) = Compound(rand(N...),zeros(N...,length(N)),ones(N...,length(N),length(N)))

julia> Compound((2,2)).vector_field
2Ă—2Ă—2 Array{Float64,3}:
[:, :, 1] =
 0.0  0.0
 0.0  0.0

[:, :, 2] =
 0.0  0.0
 0.0  0.0
2 Likes
13

I’m glad you found a solution that works for you! I do have some more thoughts reading through all the replies though.

What aspect is ugly? The definition itself, the REPL printout Compound{5, 6, 7}, or something else? If it is the definition, then of course you’re out of luck with this method. If it’s the printout or having to refer to a type with 3 parameters elsewhere in the code, then there are simple enough remedies.

If the N is necessary for the method, you would still need to write it this way I think. If not, you can omit it either way. Note also that you can write parameterized types with only some of the parameters defined. The parameters are interpreted in order.

julia> AbstractCompound == AbstractCompound{N} where N
true

julia> Compound{1}
Compound{1,M,P} where P where M

julia> Array{Int} # same idea as this
Array{Int64,N} where N

If for whatever reason you need the type to be totally concrete for some usage but don’t want to refer to it with all 3 parameters, you could define something like compoundtype(n) = Compound{n, n+1, n+2} as sugar for referring to the appropriate type. You could also define aliases like const Compound1 = Compound{1,2,3} which is similar to what you’re doing now.

Finally, the constructor you defined for Compound(N) would work (as written) just as well for a type with 3 parameters under the hood.

All of this is not to say that metaprogramming and having multiple types for Compound1, Compound2, etc. is not also a perfectly good way to go!

2 Likes
14

Nice! I didn’t realize you could drop some of the parameters. That does clean up the usage a good deal.

I still prefer strictly defining n,n+1,n+2 if possible since there really is only one parameter for the type, not three.

15

Point of fact - I ended up doing this since it was the closest to what I had already, so why reinvent the wheel?

1 Like