Reducing allocations

New to Julia
21

This blog post could be useful.

22

But why not pass this argument into the function as an SMatrix from the beginning. Everything should be an SArray from the beginning.

23

atoms.lVecs needs to be an MMatrix because the values inside will change. Are you saying that I should convert the MMatrix to an SMatrix before passing it in?

Side question: I notice that converting to an SMatrix /MMatrix can be done in several different ways, some allocate and some do not. What is the preferred way to do it (@SMatrix, SMatrix{3,3,Float64,9}(atoms.lVecs))

24

Not necessarily. Are you really sure they need to change? Can’t you just replace the matrix instead of updating it? I use StaticArrays.jl a lot, but I never use MArrays, because I can always just replace instead of update.

Give it an extra go, maybe you will come to the conclusion that mutation isn’t necessary after all.

3 Likes
25

I love that idea.. a couple of questions if you will. Below is my struct that stores the lVecs

struct atoms
    title::String
    latpar::Float64
    lVecs::MMatrix{3,3,Float64,9}
    coordSys::Vector{String} # 'D' for Direct or 'C' for cartesian
    positions:: Vector{SVector{3,Float64}}  # List of all the atomic basis vector separated by type 
    velocities:: Vector{SVector{3,Float64}} # Velocities of the atoms
    masses:: Vector{Float64}
    atomTypes:: Vector{Int64} # Integers representing the type for each basis atom
    species:: Vector{String}  # Species of the atoms 
    energies::MVector{2,Float64}
    lj_vec::Vector{Float64} # Vector of 1/r^6 and 1/r^12 values for each unique pair of atom types
end

Since this isn’t a mutable struct I won’t be able to replace lVecs with a new version. I think the speed associated with working with an immutable struct (albeit with MVectors inside) is greater than working with a mutable struct even though that would allow lVecs to be an SVector. Am I wrong?

26

Can you have two structs, mutable and immutable?

27

I would like to know if your new struct that is not mutable but does have a String, and a couple of Vectors of Strings as well as an unsized Vector of Float64’s will live on the heap. My guess is that the compiler does not know the size of these so it will be on the heap, but I would like someone with much more expertise than I to comment on this.

28

Me too, good question. Just some context: coordSys really is just a single char, but it will change between two different chars, that’s why i put it into a vector. It felt like a hack but I didn’t see any other way. species won’t ever need changing, but can you make an SVector{2,String}?

29

I would try replacing the entire struct, with everything in it. The Vectors can be passed by reference, and the StaticArrays can be replaced.

30

That feels awkward to me especially since I would be copying most of the struct fields by reference and could easily make an unintended modification. And I ran a little experiment to compare repeated modifcation of an MMatrix vs repeated creation of an SMatrix:

using StaticArrays
using BenchmarkTools

# Function to modify MMatrix in-place
function modify_mmatrix!(M::MMatrix{3,3,Float64})
    for _ in 1:1000
        for i in 1:3, j in 1:3
            M[i,j] = rand()  # Simulate modification
        end
    end
end

# Function to create new SMatrix each iteration
function modify_smatrix(S::SMatrix{3,3,Float64})
    for _ in 1:1000
        S = SMatrix{3,3}(rand(3,3))  # Simulate replacement
    end
   # return S
end

# Initialize matrices
M = MMatrix{3,3}(rand(3,3))
S = SMatrix{3,3}(rand(3,3))

# Benchmark
println("MMatrix in-place modification:")
@btime modify_mmatrix!($M) #(25 micro seconds, 0 allocations)

println("SMatrix replacement:")
@btime modify_smatrix($S). # (67 microseconds, 2000 allocations)

It seems that all those modifications of the MMatrix are “in-place” and hence done on the stack. If I ever want to do matrix multiplication or something I can just convert my MMatrix to an SMatrix first. What say ye?

31

It’s definitely true that SMatrix is friendlier to the compiler optimizer than MMatrix (for small-enough size, like here). It’s more transparent to the compiler, making it more likely to get optimized well.

32

But Mmatrix performs better in this example.

33

That’s not how to create a random SMatrix. In your example you create a regular 3x3 Array{Float64, 2} in each iteration, and afterwards you convert it.

Instead, directly create an SMatrix:

S = rand(SMatrix{3,3,Float64,9})

(Not sure if all the parameters are needed for performance, test that).

1 Like
34

Even so, this particular benchmark might not be ideal, since I suspect that the rand part might be dominating the runtime, rather than the object creation/update. Also try a more lightweight operation like adding or multiplying with a number.

35

Do you the SetField.jl might be useful for your code?

36

I believe you mean Accessors.jl

37

My instinct would be to avoid MMatrix because using them often doesn’t lead to the desired result (no allocations) and their behavior is not easy to reason about. On other hand, I don’t think it’s too important whether Atoms is allocated in heap or not as it seems to be a singular object containing large vectors (that will be in the heap in any case). I’d imagine making it mutable should not change much? If you’d rather keep it immutable, one option is to have lVecs::Ref{SMatrix{3,3,Float64,9}}, essentially a vector of one.

Fun fact, even the bracket syntax produces an SMatrix with SVector inputs: new_cell = [new_vec l2 l3].

38

Your experiment does not seem meaningful:

  • It’s comparing apples to oranges, as modify_smatrix and modify_mmatrix! do different things. modify_mmatrix! mutates the argument, while modify_smatrix is a no-op. I’m not talking about the implementation here, just about the semantics.

  • As pointed out by DNF, the SMatrix function is not implemented well, resulting in unnecessarily bad performance.