Parallelising a for loop

New to Julia
1

I need to run a simple for loop in parallel; however, my attempts at using Distributed for loops and Threads have not worked. The code is attached below and I wish to parallelise the inner for loop, it is noted in the comments.

using LinearAlgebra, BenchmarkTools

function multisensorKalmanFilter()
# Admin
xDimension = 6
zDimension = 3
n = 1000
simulationLength = 3000
# Motion model
T = 20e-3
A = [Array(1.0I, 3, 3) T*Array(1.0I, 3, 3); zeros(3, 3) Array(1.0I, 3, 3)]
ATranspose = A'
u = zeros(xDimension)
u[3] = -0.5*(9.81)*(T^2)
u[6] = -(9.81)*(T)
# Process noise
R = 4*I
# Measurement models
numberOfSensors = 8
C = [Array(1.0I, 3, 3) zeros(3, 3)]
CTranspose = C'
Q = [I*(i/4) for i = 1:numberOfSensors]
# Groundtruth
groundTruth = zeros(xDimension)
groundTruth[1:2] = [-100; -100] .+ 2*[100; 100].*rand(2)
groundTruth[4:5] = [10; 10] .+ 4*randn(2)
groundTruth[6] = 20 .+ 2*randn()
# State variables
mu = groundTruth + randn(xDimension)
S = 4*Array(1.0I, xDimension, xDimension)
# Output
meanEstimate = repeat(mu, outer=(1, simulationLength))
# Recursive state estimation
for i = 1:simulationLength
    # Predict the targets state
    groundTruth = A*groundTruth + u
    muPredicted = A*mu + u
    SPredicted = A*S*ATranspose + R
    # Measurment updates
    hGlobal = zeros(xDimension, numberOfSensors)
    KGlobal = zeros(xDimension, xDimension, numberOfSensors)
    # This for loop can be run in parallel
    for j = 1:numberOfSensors
        # Determine local posterior
        measurement = C*groundTruth + Q[j]*rand(zDimension)
        G = SPredicted*CTranspose/(C*SPredicted*CTranspose + Q[j])
        muPosteriorLocal = muPredicted + G*(measurement - C*muPredicted)
        KPosteriorLocal = inv((I - G*C)*SPredicted)
        # Convert to canonical form
        KGlobal[:, :, j] = KPosteriorLocal
        hGlobal[:, j] = KPosteriorLocal*muPosteriorLocal
    end
    # Merge the posterior distributions
    KPredicted = (1 - numberOfSensors)*inv(SPredicted)
    hPredicted = KPredicted*muPredicted
    KMerged = KPredicted + reshape(sum(KGlobal, dims=3), (xDimension, xDimension))
    hMerged = reshape(sum(hGlobal, dims=2), xDimension)
    # Posterior estimate
    S = inv(KPredicted +  KMerged)
    mu = S*(hPredicted + hMerged)
    meanEstimate[:, i] = mu
end
return meanEstimate
end

@btime meanEstimate = multisensorKalmanFilter()

The code implements a multisensor parallel update Kalman filter, my hope is to get it to run as fast as possible. I’m new to Julia, are Julia’s parallel and multi-threading capabilities appropriate for such a task?

2

The inner loop could potentially see some speedup with threading. Threading does not make a huge difference for code that allocates a lot of memory. Try making use of the insight gained from your other question and see if you can apply it here (I think you can).

BTW, Q*rand should probably call randn?

1 Like