CODEWITHSUNDEEP
multithreadingmatlabparallel-processingjuliafinite-element-analysis
nicolò
nicolò

Reputation: 31

Is the linear system solver \ also multi threaded in Julia as in Matlab? And how to "multithread" it in Julia?

I am trying to compare speed and performance between Matlab and Julia. I am looking at a code that does topology optimization of a continuum structure subjected to a given load. The code I am looking at is the public code topopt88.m: https://www.topopt.mek.dtu.dk/Apps-and-software/Efficient-topology-optimization-in-MATLAB

Essentially it is an iterative algorithm where in every iteration a Ax=b system is solved (x=A\b), where A depends on the structural design (it is the finite element stiffness matrix) and it is updated in every iteration.

In Julia the same code runs slower than Matlab. I have done some code optimization in Julia, declaring types in function definitions, using functions as much as possible, avoiding global variables, and implementing other tips I found in the internet. But Julia is still slower than the same Matlab code (same in the sense of conceptual steps).

My question: since Matlab system solve "\" is multi threaded by default, is it true the same for Julia? If not, how to multi thread Julia's \ operator, or to get speed-ups from parallelization similarly?

Upvotes: 3

Views: 858

Answers (2)

nicolò
nicolò

Reputation: 31

The following is the current best solution that I found and that gives the best performance to me. Outside Julia in a terminal in Windows type set JULIA_NUM_THREADS=4. Start Julia and use the Linear Algebra module this way: using LinearAlgebra BLAS.set_num_threads(1)

In this way there will be 4 threads, and by setting the BLAS threads to 1 it basically will not use its own pool that seems to conflict with Julia threads pool, or at least to slow the performance for me. I got some hint in this regard here: https://discourse.julialang.org/t/julia-threads-vs-blas-threads/8914/15 My current BLAS setting in Julia 1.5.1 is julia> BLAS.vendor() :openblas64

Using Pardiso.jl with the default MKL sparse solver didn't help. I still need to figure out if there is a way to increase performance there too and how.

Upvotes: 0

Przemyslaw Szufel
Przemyslaw Szufel

Reputation: 42244

By default Julia is using BLAS/OpenBLAS and you can configure it to be multi-threaded. This requires running Julia in a multi-threaded setting and setting the BLAS.set_num_threads() configuration.

Here is how:

Before starting Julia:

set JULIA_NUM_THREADS=4

or on Linux

export JULIA_NUM_THREADS=4

Now lets test in a single thread:

julia> using BenchmarkTools, Random

julia> const b = rand(50);

julia> const A = rand(50,50);

julia> @btime A \ b;
  424.899 μs (4 allocations: 20.61 KiB)

Now the multi-threaded:

julia> using LinearAlgebra

julia> BLAS.set_num_threads(4)

julia> @btime A \ b;
  175.701 μs (4 allocations: 20.61 KiB)

You can see that on my machine the speedup is over 2x.

Yet another option for speedup would be to switch Julia to MKL.

julia> pkg"add https://github.com/JuliaComputing/MKL.jl"
julia> pkg "build MKL"
# restart Julia

I have not been using MKL.jl so if you give it a try, please write how it benchmarked.

Upvotes: 3

Related Questions