Immutable FieldVectors from StaticArrays.jl

[ranocha]

SArray from StaticArrays works fine, but using a custom immutable type which is a subtype of FieldVector, recursivecopy! fails.

using StaticArrays
using RecursiveArrayTools

a = [SVector(1.,2.)]
b = similar(a)
recursivecopy!(b, a)


struct Bad <: FieldVector{2,Float64}
    a::Float64
    b::Float64
end

a = [Bad(1.,2.)]
b = similar(a)
recursivecopy!(b, a)

This results in

type Bad is immutable

Stacktrace:
 [1] copy!(::Bad, ::Bad) at ./multidimensional.jl:805
 [2] recursivecopy!(::Array{Bad,1}, ::Array{Bad,1}) at /.../.julia/v0.6/RecursiveArrayTools/src/utils.jl:21
 [3] include_string(::String, ::String) at ./loading.jl:515

According to the documentation of StaticArrays.jl

it is also worth noting that FieldVectors may be mutable or immutable

Thus, the following works fine.

mutable struct NotBad <: FieldVector{2,Float64}
    a::Float64
    b::Float64
end
a = [NotBad(1.,2.)]
b = zeros(a)
recursivecopy!(b, a)

I came across this while I was using OrdinaryDiffEq.jl with u of the form [Bad(1.,2.)]. However, this has been working some months ago, see SciML/OrdinaryDiffEq.jl#50.

[ChrisRackauckas]

I see, this needs to go down the path for Number.

[ChrisRackauckas]

Is it not taking this path?

https://github.com/JuliaDiffEq/RecursiveArrayTools.jl/blob/master/src/utils.jl#L9

[ranocha]

No, @whichtells me that it takes the path https://github.com/JuliaDiffEq/RecursiveArrayTools.jl/blob/master/src/utils.jl#L19

[ChrisRackauckas]

I just changed master to use StaticArray instead of SArray. Trying it now.

[ChrisRackauckas]

Yeah, that works, but I think then it's broken for mutable FieldVectors. I am not sure how to query for mutability.

[ranocha]

Yeah, it works fine for me on master. Thank you!

[ChrisRackauckas]

Does mutable FieldVector work? Or will they alias after that?

[ranocha]

They alias.

a = [NotBad(1.,2.)]
b = zeros(a)
recursivecopy!(b, a)

@show a
@show b

a[1][1] = 5.
@show a
@show b

yields

a = NotBad[[1.0, 2.0]]
1-element Array{NotBad,1}:
 [1.0, 2.0]
b = NotBad[[1.0, 2.0]]
1-element Array{NotBad,1}:
 [1.0, 2.0]

a = NotBad[[5.0, 2.0]]
1-element Array{NotBad,1}:
 [5.0, 2.0]
b = NotBad[[5.0, 2.0]]
1-element Array{NotBad,1}:
 [5.0, 2.0]

[ChrisRackauckas]

@andyferris do you know of a way to handle this correctly? The problem is one needs to use copy! while the other needs to a[i] = b[i], but I am not sure how to, at compile time, check for this.

[andyferris]

Yeah - we need a trait (in Base) for array mutability (answering the question: "does setindex! work?).

I'm happy to support a PR to Base, or to StaticArrays, or for a 3rd package to include the trait and include it in StaticArrays. Otherwise, I'm sure I will get around to it (one day ;) )

I've generally just assumed that static arrays are immutable since creating a copy is typically cheap(ish).

[ChrisRackauckas]

I've generally just assumed that static arrays are immutable since creating a copy is typically cheap(ish).

I guess copy of a mutable FieldVector is a good short term fix.

I'm happy to support a PR to Base, or to StaticArrays, or for a 3rd package to include the trait and include it in StaticArrays. Otherwise, I'm sure I will get around to it (one day ;) )

Yeah, we should really make sure this is at least part of 1.0.

[andyferris]

Yeah, we should really make sure this is at least part of 1.0.

100% agree.

[ranocha]

Could the function is_immutable_typefrom https://stackoverflow.com/questions/40706962/how-to-check-a-type-for-mutability solve the problem?

using StaticArrays

struct ImmutableFV <: FieldVector{2,Float64}
    a::Float64
    b::Float64
end

mutable struct MutableFV <: FieldVector{2,Float64}
    a::Float64
    b::Float64
end

Base.@pure is_mutable_type(x::DataType) = x.mutable

@show is_mutable_type(ImmutableFV)
@show is_mutable_type(MutableFV)

# small check whether is_immutable_type is determined at compile time
function two_branches(x)
    if is_mutable_type(typeof(x))
        1
    else
        1.
    end
end

immutable_fv = ImmutableFV(1.,2.)
mutable_fv = MutableFV(1.,2.)

@code_warntype two_branches(immutable_fv)
@code_warntype two_branches(mutable_fv)

This returns

julia> @code_warntype two_branches(immutable_fv)

Variables:
  #self#::#two_branches
  x::ImmutableFV

Body:
  begin 
      goto 3 # line 4:
      3:  # line 6:
      return 1.0
  end::Float64

julia> @code_warntype two_branches(mutable_fv)
Variables:
  #self#::#two_branches
  x::Any

Body:
  begin  # line 4:
      return 1 # line 6:
  end::Int64

[andyferris]

The primary concern is getting this to correctly reflect whether setindex! will work.

For instance, RowVector is an immutable wrapper but will setindex! if the parent array can. The trait needs to reflect this.

[ranocha]

Ah, okay, I didn't think of something like this.

For StaticArrays, setindex! uses setfield!, which is a built-in function. method_exists does not seem to work with setfield!, it just returns true.