Pointer-like field in a struct?

General Usage
1

Hi,
I hope somebody here can help me with this one.
I’m building a Monte Carlo code for a system of particles. Position coordinates of these N particles are stored in a 3xN array of coordinates r[:,:], like for instance for N=5
‘’’
0.0719552 0.820621 0.169617 0.336838 0.620763
0.218954 0.198641 0.793058 0.661784 0.32143
0.448448 0.377976 0.533512 0.278445 0.821139
‘’’
where each particle sits in a different column, and rows indicate x,y and z for that particle, respectively.

Now I need to define an object (called a bead) for each of these particles, which I want to define in a mutable struct.That’s because every particle has some simulation parameters corresponding to it, and its own coordinates, together with the coordinates of its left particle, and the coordinates of its right particle. Then, along the simulation, the particle positions change according to a Metropolis criteria. What I would like to do, however, is to perform these coordinates changes in r[:,:], but such that the changes automatically apply to the coordinates in the beads. I have tried to define the beads as in
‘’’
mutable struct Bead{T}
δt :: T;
rleft :: Array{T};
rright :: Array{T};
Bead{T}() where T = new();
end
‘’’
so that I can do
‘’’
kind=Float32
b1 = Bead{kind}();
b2 = Bead{kind}();
b1.δt = b2.δt = 0.01;
b1.rleft = view(r,:,2);
b1.rright = view(r,:,3);
b2.rleft = view(r,:,3);
b2.rright = view(r,:,4);
‘’’
hoping that, being a view what I use, changes in r[:,:] are automaticaly updated to b1.rleft, b1.rright, b2.rleft and b2.rright.
But that does not really work, and I think that this is becasue I define
rleft and rright in the Bead struct as arrays, which copy the value they get
at the point I initialize them, but do not get updated when r[:,:] is updated.
If I was writting the code in C, I could define the rleft and rright fields in the Bead struct as pointers to r[:,:], which would then solve the problem. But in Julia… I don’t know.

So the question is: is there a way to implement this pointer-like behaviour?
Something that would automatically update the rleft and rright fields in b1 and b2 when changes in r[:,:] are prodced…

Thanks in advance,

Ferran.

2

This should already work. Julia uses pass by sharing rather than pass by value. Values do not get copied unless you copy them.

3

Uhmm… don’t get it to work.
Here’s what I tried:

mutable struct Bead{T}
    δt     :: T;
    rleft  :: Array{T};
    rright :: Array{T};
    Bead{T}() where T = new();
end

kind = Float32;

# define an array of positions from part of a chain
# 
r  = rand(3,5)

3×5 Matrix{Float64}:
 0.74328   0.253234  0.506628  0.581835  0.0087231
 0.612277  0.934715  0.723907  0.380912  0.537104
 0.374414  0.741888  0.985717  0.174058  0.0596843

b1 = Bead{kind}();
b2 = Bead{kind}();

b1.rleft  = view(r,:,2);
b1.rright = view(r,:,3);
b2.rleft  = view(r,:,3);
b2.rright = view(r,:,4);

hcat(b1.rright, b2.rleft)

3×2 Matrix{Float32}:
 0.506628  0.506628
 0.723907  0.723907
 0.985717  0.985717

r[:,3] = ones(3)
display(r)

3×5 Matrix{Float64}:
 0.74328   0.253234  1.0  0.581835  0.0087231
 0.612277  0.934715  1.0  0.380912  0.537104
 0.374414  0.741888  1.0  0.174058  0.0596843

hcat(b1.rright, b2.rleft)

3×2 Matrix{Float32}:
 0.506628  0.506628
 0.723907  0.723907
 0.985717  0.985717

so you see, the coordinates in b1.rleft etc are not updated to the new values :frowning:

Best,

Ferran.

​

4

If you change Bead to

mutable struct Bead{T}
    δt     :: T;
    rleft  :: AbstractArray;
    rright :: AbstractArray;
    Bead{T}() where T = new();
end

it would work. The problem is that there is no way to make a subarray of a Matrix{Float64} fit into an `Array{Float32}

1 Like
5

Why not make the abstract array type a parameter?

julia> mutable struct Bead{T, A <: AbstractArray{T}}
           δt     :: T
           rleft  :: A
           rright :: A
           function Bead(rleft::A, rright::A) where {T, A <: AbstractArray{T}}
               new{T,A}(zero(T), rleft, rright)
           end
       end

julia> kind = Float32;

julia> r = rand(kind, 3, 5);

julia> b1 = Bead(view(r,:,2), view(r,:,3))
Bead{Float32, SubArray{Float32, 1, Matrix{Float32}, Tuple{Base.Slice{Base.OneTo{Int64}}, Int64}, true}}(0.0f0, Float32[0.29957736, 0.95468074, 0.19946498], Float32[0.6553544, 0.96422833, 0.7301229])

julia> b2 = Bead(view(r,:,3), view(r,:,4))
Bead{Float32, SubArray{Float32, 1, Matrix{Float32}, Tuple{Base.Slice{Base.OneTo{Int64}}, Int64}, true}}(0.0f0, Float32[0.6553544, 0.96422833, 0.7301229], Float32[0.1161111, 0.386643, 0.23078555])

julia> hcat(b1.rright, b2.rleft)
3×2 Matrix{Float32}:
 0.655354  0.655354
 0.964228  0.964228
 0.730123  0.730123

julia> r[:,3] .= 1
3-element view(::Matrix{Float32}, :, 3) with eltype Float32:
 1.0
 1.0
 1.0

julia> r
3×5 Matrix{Float32}:
 0.229182  0.299577  1.0  0.116111  0.27342
 0.352749  0.954681  1.0  0.386643  0.952313
 0.993329  0.199465  1.0  0.230786  0.232765

julia> hcat(b1.rright, b2.rleft)
3×2 Matrix{Float32}:
 1.0  1.0
 1.0  1.0
 1.0  1.0

To clarify, the original issue is that you have declared the field types to be an Array{T}. view creates a SubArray{T}. When assigning to a field, Julia will implicitly convert the right hand side type to the field type. In doing so, it makes a copy of the array.

To avoid making the copy, we need to make the field type compatible with the type being assigned so that no conversion and thus copying is done. Oscar’s method declares the field to be an abstract type.

In my example, I use a type parameter. This allows the field type of rleft and rright to be concretely determined from the concrete type of the Bead.

2 Likes
6

Hi guys,

Great answers! I’ll try both and see, probably tomorrow morning.
And btw, I did not now that struts could use parameters, which is something I thought about but believed to not be possible. Thanks for that piece of info, too :slight_smile:

Best regards,

Ferran

7

Note that this is a classic case where you should consider StaticArrays.jl, i.e. use a length-N array of 3-component SVectors, e.g. Vector{SVector{3,Float64}}. Under the hood, the storage is the same as a 3xN array, but SVector is typically much easier to work with, and often faster because all of the length-3 loops will be unrolled.

e.g. you can just do r[3] = r[1] + r[2] set the third vector to the sum of the first two, and it allocates no temporary arrays on the heap and all the loops are unrolled: it just does 3 additions and an assignment. In contrast, r[:,3] = r[:,1] + r[:,2] allocates a temporary array (unless you use dot operations), and also involves a runtime loop (because the compiler doesn’t know that the loop is of length 3).

And if you want to do the same operation to all of the vectors, you can use broadcasting, e.g. r .= update.(r) will call r[i] = update(r[i]) on every element of r, and since these are SVectors the update function can just use natural vector operations and doesn’t need to worry about allocations of 3-component vectors.

2 Likes
8

Hi again,

yesterday I found the time to test the proposals and thanks, yes, they work well :slight_smile:
One thing that I liked was the parametrized version proposed by @mkitti, but there is something I’d like to avoid. And that concerns the initialization at the end of the deffinition

function Bead(rleft::A, rright::A) where {T, A <: AbstractArray{T}}
               new{T,A}(zero(T), rleft, rright)

as it it, the new() part has to say something about all fields in the struct. If you look at my original (non-working) version, I just write new() because I prefer to
initialize things externally. I do that because most of my structs undergo field changes with time and I do not want to haveing to care about that every time.
Yes I know this may not be the best practica, but it has worked for me so far and, honestly, I do not want to change that…
So I was wondering if there is a way to just specify rleft and rright at initialization, or not specify anything about them until the object is created externally, or…
Maybe you know a way to do that?

Best regards and thanks,

Ferran.