# How to make a vector of parametric composite types?

So I am not sure exactly what you want. Maybe instead you are looking for something like this?

``````@enum Material Sic Al Vac

immutable Layer
material::Material
thickness::Float64
end

permittivity(l::Layer) = permittivity(Val{l.material})

function permittivity(::Type{Val{Al}})
println("The permittivity of Aluminum is grand indeed")
end

function permittivity(::Type{Val{Sic}})
println("But even grander is Silicium")
end

#Then
l = Layer(Al, 4)
permittivity(l)
l2 = Layer(Sic, 4)
permittivity(l2)
``````
``````@time ml = [Layer(Al,0.0) ; Layer(Vac,0.0)]
0.000008 seconds (7 allocations: 336 bytes)
``````
3 Likes

Woow, thanks a lot!
This is really nice.
It will make the code much more user friendly.

Thanks again
Olivier

1 Like

I think this will cause calculating the `permitivity` to be slow, since it will be doing dynamic (runtime) dispatch.

Would it be possible to use a dictionary for the permitivities instead?

1 Like

You’re right.

``````@enum Material Sic Al Vac

immutable Layer
material::Material
thickness::Float64
end

permittivity(l::Layer) = permittivity(Val{l.material})
permittivity (generic function with 1 method)

function permittivity(::Type{Val{Sic}})
2
end

using BenchmarkTools

@benchmark permittivity(l2)
BenchmarkTools.Trial:
memory estimate:  16 bytes
allocs estimate:  1
--------------
minimum time:     331.807 ns (0.00% GC)
median time:      335.516 ns (0.00% GC)
mean time:        346.447 ns (0.45% GC)
maximum time:     8.241 μs (93.95% GC)
--------------
samples:          10000
evals/sample:     223
time tolerance:   5.00%
memory tolerance: 1.00%
``````

vs

``````abstract OptProp

immutable Al <: OptProp end

immutable Layer{T <: OptProp}
material::Type{T}
thickness::Float64
end

function permittivity{T<:Al}(x::Layer{T})
2
end
l1 = Layer(Al, 0.)

using BenchmarkTools
@benchmark permittivity(l1)
BenchmarkTools.Trial:
memory estimate:  0 bytes
allocs estimate:  0
--------------
minimum time:     12.587 ns (0.00% GC)
median time:      12.725 ns (0.00% GC)
mean time:        12.890 ns (0.00% GC)
maximum time:     46.025 ns (0.00% GC)
--------------
samples:          10000
evals/sample:     999
time tolerance:   5.00%
memory tolerance: 1.00%
``````

Thanks! (and sorry @Olivier_Merchiers)

1 Like

@Olivier_Merchiers can you give a little more info - e.g. why do you use dispatch on the permittivity function?

Thanks for this correction.

But there remains the issue with allocation and time for the vector of layers with this implementation:

``````abstract OptProp

immutable Al <: OptProp end
immutable Vac <: OptProp end

immutable Layer{T <: OptProp}
material::Type{T}
thickness::Float64
end

@benchmark ml = [Layer(Al,0.0); Layer(Vac,0.0)]
BenchmarkTools.Trial:
memory estimate:  1.97 kb
allocs estimate:  48
--------------
minimum time:     10.937 μs (0.00% GC)
median time:      11.656 μs (0.00% GC)
mean time:        11.802 μs (0.00% GC)
maximum time:     72.133 μs (0.00% GC)
--------------
samples:          10000
evals/sample:     1
time tolerance:   5.00%
memory tolerance: 1.00%

``````

While with

``````abstract OptProp

immutable Al <: OptProp end
immutable Vac <: OptProp end

immutable Layer
material::OptProp
thickness::Float64
end

@benchmark ml = [Layer(Al(),0.0); Layer(Vac(),0.0)]
BenchmarkTools.Trial:
memory estimate:  160.00 bytes
allocs estimate:  3
--------------
minimum time:     43.866 ns (0.00% GC)
median time:      47.852 ns (0.00% GC)
mean time:        61.659 ns (21.21% GC)
maximum time:     2.199 μs (96.82% GC)
--------------
samples:          10000
evals/sample:     990
time tolerance:   5.00%
memory tolerance: 1.00%

``````

This is quite important for me, since I need to pass such multilayer structures to other functions that will compute reflectivity and other physical quantities.

Finally, doing the test with the permittivity function

``````function permittivity(x::Al)
2
end

``````

I get

``````@benchmark permittivity(Al())
BenchmarkTools.Trial:
memory estimate:  0.00 bytes
allocs estimate:  0
--------------
minimum time:     1.409 ns (0.00% GC)
median time:      1.680 ns (0.00% GC)
mean time:        1.711 ns (0.00% GC)
maximum time:     90.449 ns (0.00% GC)
--------------
samples:          10000
evals/sample:     1000
time tolerance:   5.00%
memory tolerance: 1.00%

``````

Which is almost 10 times faster.

@mkborregaard Sorry, I missed this last question before submitting the other reply.

Yes of course.
I need to make reflectivity calculations for many different materials. These materials have optical properties that are defined by their permittivity function.
With every material corresponds a unique permittivity function which returns a complex number.
This number is used to compute the reflectivity.

To keep my code generic, and not using if-then-else in a big permittivity function containing all possible definitions for each material, I wanted to use multiple dispatch to get the right permittivity method corresponding with each material.

To give some examples:
For vacuum this would be

``````permittivity(material::Vacuum) = 1.0 +im*0.0
``````

For frequency dependent (hence w in the parameter list) optical properties you would have something like

``````function permittivity(material::Cbn,w)
eps_fin = 4.46 + 0.0*im
return eps_fin*(w^2-w_lo^2 + im*gamma*w)/(w^2-w_to^2 + im*gamma*w)
end
``````

Let me know if this explanation is not clear

1 Like

That’s why I suggested to use a dictionary instead of the big `if-else`.

How would you implement such a Dict?

If the permittivity for each material was a constant value, I think I could build such a dict.
But this is not the case here, since the permittivity for each material is a function of the frequency as shown in this example:

``````function permittivity(material::Cbn,w)
eps_fin = 4.46 + 0.0*im
return eps_fin*(w^2-w_lo^2 + im*gamma*w)/(w^2-w_to^2 + im*gamma*w)
end
``````

where w is the frequency in rad/s and can take any real value between 0 and infinity.

That would imply a dict of functions?

Well maybe your original model is best, I really don’t see what is wrong with it But the reason your last permittivity function is fast is because there is no dispatch at all.
You could have concrete parametric functions in your first example like

``````    abstract OptProp

immutable Al <: OptProp end
immutable Vac <: OptProp end

immutable Layer{T <: OptProp}
material::T
thickness::Float64
end

# and then
permittivity{T <: Al}(l::Layer{T}) = 2

# or the clearer
permittivity(l::Layer) = permittivity(l.material)
permittivity(x::Al) = 2
``````

Note that in benchmark it is important how to set it up for comparison

``````l = Layer(Al(), 0.)
@benchmark permittivity(l)
BenchmarkTools.Trial:
memory estimate:  0 bytes
allocs estimate:  0
--------------
minimum time:     15.057 ns (0.00% GC)
median time:      15.652 ns (0.00% GC)
mean time:        15.877 ns (0.00% GC)
maximum time:     46.872 ns (0.00% GC)
--------------
samples:          10000
evals/sample:     998
time tolerance:   5.00%
memory tolerance: 1.00%

@benchmark permittivity(Layer(Al(), 0.))
BenchmarkTools.Trial:
memory estimate:  0 bytes
allocs estimate:  0
--------------
minimum time:     1.500 ns (0.00% GC)
median time:      1.502 ns (0.00% GC)
mean time:        1.506 ns (0.00% GC)
maximum time:     3.348 ns (0.00% GC)
--------------
samples:          10000
evals/sample:     1000
time tolerance:   5.00%
memory tolerance: 1.00%
``````
2 Likes

You need to interpolate:

``````@benchmark permittivity(\$l)
``````

This then gives the same time.

1 Like

See here for details about interpolation.

1 Like

This seems to be the best solution.

Thanks, this is great - in trying to help I also (re-)learn some useful things myself. I only hope this hasn’t created too much noise for anyone reading this thread

1 Like

Thanks a lot for your help @mkborregaard, @dpsanders.

This was really instructive. Things are much clearer now to me.

With the interpolation, the benchmark results for the vector of Layers become also identical for the different approaches even with the enum type approach. Too bad the enum types are slow for the permittivity calculation. I really liked them

1 Like

If it was the syntax you liked about the enum type approach you could always define a little macro to define the abstract and immutable types:

``````macro myenum(args...)
exs = Expr[]
push!(exs, :(abstract \$(args[1])))
for arg in args[2:end]
push!(exs, :(immutable \$(arg) <: \$(args[1]) end))
end
esc(Expr(:block, exs...))
end
``````

The enum example looked like this:

``````@enum Material Sic Al Vac

immutable Layer
material::Material
thickness::Float64
end

permittivity(l::Layer) = permittivity(Val{l.material})

function permittivity(::Type{Val{Sic}})
println("pompuous nonsense")
end
``````

With the myenum macro you could write the last, recommended, approach as

``````@myenum Material Sic Al Vac

immutable Layer{T <: Material}
material::T
thickness::Float64
end

permittivity(l::Layer) = permittivity(l.material)

function permittivity(x::Al)
println("pompuous nonsense")
end
``````

Which I actually like even better.

1 Like

If I understand correctly, I think you have settled with solution:

``````abstract OptProp

immutable Al <: OptProp end
immutable Vac <: OptProp end

immutable Layer
material::OptProp
thickness::Float64
end
``````

I think this is a good idea. I don’t see a speed issue here (I will elaborate below).

I do suggest you change your names, though. I noticed @mkborregaard thought you were talking about aluminum & silicium, but my guess is that you were talking about alumina and silicon carbide (SiC) - otherwise you would have just called your type Si (unless you were worried of confusions with the french “if” ).

I would also rename `OptProp` to be `Material`… since this `Al` & `Vac` are not optical properties, but materials.

1 Like

### `Material{:Alumina}` solution

(vs: `immutable Alumina < Material; end`)

In your particular case, I don’t believe there is anything significantly different (performance-wise) with the Julia implementation the `Material{}` solution than with the `immutable Alumina` solution.

I use both solutions quite often myself. I often prefer the `immutable Alumina` solution simply because it avoids mistakes:

You can accidentally construct objects with typos using parametric types, unless you are careful to raise exceptions when the wrong one is created:

``````m = Material{:Alumina}() #Great
m = Material{:Aluminq}() #You get a silent error here
``````

This is not an issue using the `immutable Alumina` solution:

``````m = Alumina() #Great
m = Aluminq() #Julia will raise an error (very helpful)
``````

### Type stability & layers

When you define your `Layer` object as:

``````immutable Layer
material::Material
thickness::Float64
end
``````

You can now build material stacks of `Vector{Layer}` (which are type stable):

``````interfaceAB = [LayerA ; LayerB]
typeof(interfaceAB)
--> Array{Layer, 1}
``````

Now, as soon as you parameterize `Layer`, you get type instability:

``````immutable Layer{T}
material::Material{T}
thickness::Float64
end

interfaceAB = [LayerA ; LayerB]
typeof(interfaceAB)
--> Array{Layer, 1}
``````

Note that this time, `Layer` represents a “not fully concrete type” (almost like an abstract type). A type is only concrete when it is qualified with its parameters (ex: `Layer{:Alumina}`):

``````layerAl = Layer(Material{:Alumina}, 1)
interfaceAB = [layerAl; layerAl]
typeof(interfaceAB)
--> Array{Layer{:Alumina}, 1}
``````

This time `interfaceAB` will generate “type stable” code - because it is a vector of truly concrete types.

There in lies the problem: you will only get “type stable” code if all your types are concrete (ie your materials are all the same).

1 Like

### Speed Issues

Okay, now after saying all of this, what about speed?

As I just mentionned, you can only get truly “type stable” code using this methodology if you are building a stackup of the same material (sort of a useless case).

So what you are really tring to do is trying to avoid writing your own `if` statements - by leveraging Julia’s type system. I do this all the time, and I think it is great (I use this instead of using `Dict`s myself).

I will try to illustrate what changes between the 2 solutions:

``````#Parameterized:
immutable MaterialP{T}; end
typealias AluminaP MaterialP{:Alumina}
typealias VacuumP MaterialP{:Vacuum}
immutable LayerP{T}
material::MaterialP{T}
end
LayerAlP = LayerP(AluminaP())
LayerVacP = LayerP(VacuumP())

#Non-parameterized:
abstract Material
immutable Alumina <: Material; end
immutable Vacuum <: Material; end
immutable Layer
material::Material
end
LayerAl = Layer(Alumina())
LayerVac = Layer(Vacuum())
``````
``````function myfun(v::Vector{LayerP})
return [dosomething(l) for l in v]
end
function myfun(v::Vector{Layer})
return [dosomething(l) for l in v]
end
``````

Now, the following might be considered type stable, but there is an implicit “if” stmt called just before calling `dosomething()`… that is unless Julia can statically figure out what types you have (ex: if it inlines `myfun()`):

``````#Case 1
myfun([LayerAl, LayerVac])
``````

As far as I know, Julia’s type system would generate very similar code/speed if we used an array of parameterized Layers:

``````#Case 2
myfun([LayerAlP, LayerVacP])
``````

Where Julia can make its optimization, though, is if it understands that all the elements are of the same type. That way, it does not have to do its implicit `if` before calling `dosomething()`:

``````#Case 3
myfun([LayerAlP, LayerAlP])
``````

… Not often a very useful stack, though (you would probably merge these materials yourself into a single thickness).

Now here is where you loose speed with the non-parameterized system: When you create a stack `Vector` using the the same `Layer` object twice, you don’t actually get a speedup wrt Case 1:

``````#Case 4
myfun([LayerAl, LayerAl])
``````

… Again, that’s not really a useful stack, so I don’t see the point in optimizing for this.

Hope this helps.