Hello everyone,
I want to do n-body simulations with Julia. The number of particles is
small (10 to 100) but the simulation runs for tens of millions of orbits
and I want to optimize them. The most natural way to implement this is to
have one array for every particle property:
x = Array(Float64,nbodies) # x coordinates
y = Array(Float64,nbodies) # y coordinates
z = Array(Float64,nbodies) # z coordinates
u = Array(Float64,nbodies) # u coordinates
v = Array(Float64,nbodies) # v coordinates
w = Array(Float64,nbodies) # w coordinates
m = Array(Float64,nbodies) # m coordinates
That means that when I compute the force between two particles I need to
access 7 x 2 = 14 values that are not sequential in memory. So I was
thinking that it might be worth making a new type:
type TBody
x::TFloat
y::TFloat
z::TFloat
u::TFloat
v::TFloat
w::TFloat
m::TFloat
end
bodies = Array(TBody,nbodies)
So, when I compute the force between two particles, I should only need to
access two bits of memory. But on the other hand, maybe this is all a waste
of time. With only 100 particles, and 7 doubles per particle, the entire
memory cost is about 5 kb and it all fits easily inside the CPU's L1 cache.
So maybe this is just a classical example of premature optimization that
has no value whatsoever.
Can anyone give me an opinion here? Is there something I haven't thought of?
Cheers,
Daniel.
