Differential Equations.jl ODEProblem exceedingly slow

Specific Domains Numerics
1

I am trying to numerically integrate a set of master equations with DifferentialEquations.jl

I have defined the differential equation like:

function voter_model_ame!(dG,G,t,p)
    N,U  = size(G)
    for n_prime=1:N,u_prime=1:U
        u  = u_prime-1
        n = n_prime-1
        dG[n_prime,u_prime] = 0

        if n == 0
            continue
        end

        if u > n 
            continue
        end

        (u < n) && (dG[n_prime,u_prime] -= G[n_prime,u_prime]*(n-u)*(u/n))#upflip outflux
        (u > 0) && (dG[n_prime,u_prime] -= G[n_prime,u_prime]*(u)*((n-u)/n))#downflip outflux
        (u < n) && (dG[n_prime,u_prime]+= G[n_prime,u_prime+1]*((u+1)*((n-u-1)/n)))#inwards downflip flux
        (u > 0) && (dG[n_prime,u_prime]+= G[n_prime,u_prime-1]*((n-u+1)*((u-1)/n)))#inwards upflip flux
    end
    return dG 
end

This is a system of N^2 differential equations.

I can set the parameter N = 500 and I try to numerically integrate like this:

t_max = 10
p = [0]
u₀ = SAME.initialize_u0(N_individuals,p_dist)#creates initial conditions, Matrix{BigFloat} of size NxN
t_span = (0.,t_max)
prob = ODEProblem(voter_model_ame!,u₀,t_span,p)
sol = solve(prob,Tsit5(),saveat=1)

The code take ~3 seconds for N = 50 and upwards of 10min to run with N = 500.
Here are the @time return for N = 50

4.603619 seconds (77.25 M allocations: 3.787 GiB, 8.55% gc time, 46.87% compilation time)

Now this seems to be very slow and takes up an insane amount of memory.

I used the in-place style function instead of allocating a new array every time step as per the optimization instructions in the DifferentialEquations.jl I can’t use static arrays because the size of my arrays is too large.

I have a feeling I am doing something horrendously wrong but I don’t know what. There must be an efficient way to solve large systems of differential equations. What am I missing?

2

Can you @time or benchmark your rhs! function?

3 
du₀ =zeros(size(u₀))#50x50 Matrix{BigFloat}
@time SAME.voter_model_ame!(du₀,u₀,1,[0])

gives

0.001388 seconds (24.50 k allocations: 1.215 MiB)

For 10 times steps do we expect 4 orders of magnitude slowdown? It seems like it should scale linearly in the number of time steps.

Secondly, it seems like it is taking a huge amount of memory. Does ODEProblem require the use of BigFloat, arbitrary precision floats?

4

No, and in fact it should be much slower than using fixed-precision floats. Is there a specific reason you are generating u₀ with BigFloat? Using e.g. Float64 should get rid of many of the allocations in the ODE function (ideally you want this to have no allocations at all).

5

Okay wow that makes a crazy difference.

Now when I run the !rhs function I get

0.000009 seconds (1 allocation: 64 bytes)

So the BigFloat is what was slowing it down.

Thanks for the help

6

Nice, I’m glad it helped!

I’m sure that there are some important use cases for arbitrary precision floats, but from my experience, Float64 should be the go-to in virtually all situations (or if memory is scarce perhaps Float32).