For Monte Carlo simulation with same code same algorithm, how fast is Julia compared with Fortran?

@CRquantum

how long is the portion of your code that actually does the real algorithm?

You can define a macro to effortlessly unpack a large number of values at once from a single object.

Sample code:

struct Foo{A, B, C} a::A; b::B; c::C end

"""
`@unpackall_Foo(obj)` unpacks all fields of the object `obj` of type `Foo`.
"""
macro unpackall_Foo(obj)
    names = fieldnames(Foo)
    Expr(:(=),
        Expr(:tuple, names...),
        Expr(:tuple, (:($obj.$name) for name in names)...)
    ) |> esc
end

@doc @unpackall_Foo

@unpackall_Foo(obj) unpacks all fields of the object obj of type Foo .

Input:

@macroexpand @unpackall_Foo foo

Output:

:((a, b, c) = (foo.a, foo.b, foo.c))

Input:

@unpackall_Foo Foo(1, 2.0, "three")
a, b, c

Output:

(1, 2.0, "three")

With this macro, there is no need to share a large number of global constants between multiple functions.

Just pass the object foo, which contains the values shared between multiple functions, to a function as an argument each time, and unpack it in the function.

Input:

foo = Foo(1, 2.0, "three")

function f(foo::Foo)
    @unpackall_Foo foo
    @show a b c
    return
end

f(foo)

Output:

a = 1
b = 2.0
c = "three"

In the above, we unpack only a, b, and c from foo, but you can increase the number of them as much as you want.

To increase the number of fields in the style of struct Foo{A, B, C} a::A; b::B; c::C end as much as you want, you can use ConcreteStructs.jl (my example).

Jupyter notebook: public/How to define unpacking macros.ipynb at main · genkuroki/public · GitHub

3 Likes

Isn’t this just @unpack from Parameters.jl?

1 Like

The difference is as follows:

@unpackall_Foo foo
@unpack a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x, y, z = foo

It’s hard to be forced to use the latter style, so it would be great if someone clever could create a package to use the former style.

Sample code

Create a composite type with a lot of fields:

using ConcreteStructs

@concrete struct Foo a; b; c; d; e; f; g; h; i; j; k; l; m; n; o; p; q; r; s; t; u; v; w; x; y; z end

names = fieldnames(Foo)
kwargs = Expr(:parameters, (Expr(:kw, name, v) for (v, name) in enumerate(names))...)
@eval Foo($kwargs) = Foo($(names...))

foo = Foo(a = "meow", m = 'm', z = 99.99)

Output:

Foo{String, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Char, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Int64, Float64}(“meow”, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, ‘m’, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 99.99)

@unpackall_Foo case:

"""
`@unpackall_Foo(obj)` unpacks all fields of the object `obj` of type `Foo`.
"""
macro unpackall_Foo(obj)
    names = fieldnames(Foo)
    Expr(:(=),
        Expr(:tuple, names...),
        Expr(:tuple, (:($obj.$name) for name in names)...)
    ) |> esc
end

let
    @unpackall_Foo foo
    @show a b c d e f g h i j k l m n o p q r s t u v w x y z
end;
Output
a = "meow"
b = 2
c = 3
d = 4
e = 5
f = 6
g = 7
h = 8
i = 9
j = 10
k = 11
l = 12
m = 'm'
n = 14
o = 15
p = 16
q = 17
r = 18
s = 19
t = 20
u = 21
v = 22
w = 23
x = 24
y = 25
z = 99.99

Parameters.@unpack case:

using Parameters

let
    @unpack a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x, y, z = foo
    @show a b c d e f g h i j k l m n o p q r s t u v w x y z
end;

The output is the same as above.

Jupyter notebook: public/How to define unpacking macros Part 2.ipynb at main · genkuroki/public · GitHub

1 Like

If you define your struct with @with_kw, Parameters.jl will also make a macro for you which does @upack_all. See here.

Note: I would definitely consider this an anti-pattern. If your function needs all of those variables in scope, it’s a sign you need to re-factor into smaller functions.

6 Likes

Oh! I hadn’t noticed that part of the documentation! Thank you very much.

(What I’m having a little trouble with is about the compatibility of Parameters.@with_kw and ConcreteStructs.jl. ([ANN] ConcreteStructs.jl)

I just put in a fix for ConcreteStructs to allow it to work with Base.@kwdef. Parameters.@with_kw doesn’t work, though.

About 500 lines or so. It is prototype algorithm but it works very well and robust.
But it will take me some time to appreciate the logic of Julia.

Thank you very much for your suggestions.
I used Julia 1.6.1. I used 3 methods as you suggested,

rn = rand(1)[1] 

rn = rand()

const RNG = MersenneTwister(100)
rn = rand(RNG)

In my code it does not make difference. It looks like the bottleneck is not the random number generator. Julia 1.6.1 seems smart enough in running random number generator.

In my view, modern Fortran means the following:

  1. f90 and above. Use modules to arrange the whole program in a neat and clear way.
  2. implicit none, and use select kind to choose true real 8 or double precision.
  3. flexible arrays. Whenever possible prevent using hard coded array size. Set array size as a variable, use allocatable arrays when needed.
  4. Whenever possible, vectorize the code.
  5. Use types to organize objects.
  6. highly parallelized, especially MPI.

Thank you @Vasily_Pisarev @Sukera @genkuroki @DNF @lmiq @giordano @Henrique_Becker @ranocha @artkuo @Lime etc (Sorry for the @ ten user limit…) Your help and comments are very crucial and very helpful!

The reason my Julia version is 3 times slower than Fortran, is probably because those const Ref stuff in my code caused a lot of unnecessary memory operations. Hopefully by removing them the code can reach the speed of Fortran.

I think Julia is promising. I also hope Julia can be more intelligent in optimizing the code when compiling. and do more @inline, @inbound, @turbo, @unpack, etc… tricks automatically. Personally I do not mind Julia spending some more seconds in compiling as long as it can do more optimization automatically.

In one word, I wish Julia could perhaps improve in a way that it can be more difficult for the user to write slow code.

3 Likes

Personally, I imagine something like a “performance” or “debug” mode which could, for example, throw warnings (or even errors?) when using non-constant globals or, as an option, errors on allocations in a certain block of code (useful for high performance engineering and embedded systems).

Note, though, that being able to write naive / slow / bad code is also an important feature, in particular for interactive usage.

6 Likes

Given that Julia’s design encourages highly abstract code, it’s hard to even define what code is slow without running it.
Even the simplest case

function sum_itr(itr)
    s = 0
    for x in itr
        s += x
    end
    return s
end

is tricky: it is fast if itr has only integer values, and it has penalty due to type instability otherwise.

Such analysis is less of a problem when you have the end-user code where the types may be known at all call sites, but I can’t imagine a robust a priori performance analysis of a library code in a general case (but of course the library devs can include in the test suite the specific argument types they assume important / representative).

3 Likes

I really prefer my code to work first and performance can be continuously improved instead of re-writing it in a different language. It sounds like C/C++/Fortran are for you, where even the most verbose convoluted code is often not-too-slow

I hope we get a package that runs a function and just says, this is slow because of this etc. Basically an even more powerful Dispatch Analysis · JETTest.jl . It already is amazing, but I hope it can get even better.

5 Likes

I totally agree.
Having ''performance" mode and “Debug” mode is like the “Release” Mode and “Debug” mode of the intel Fortran in Visual Studio.
Wish Julia can throw more warning message when it finds some performance issues and tell user how to improve the code. In my point of view, being able to throw those ‘intelligent’ message can be really a “modern” feature for a modern language.
Also wish Julia could have its own official IDE which is easy to use and debug, this can make Julia distinguished more from things like Python too.

1 Like

Thank you very much for your efforts!

May I ask a stupid question, using the below structure as your said,

Uhm, why it will be faster than defining many variables each with the const + Ref trick as the `Fortran’ like way I did?

It seemed to me, from reading this thread, that using const + ref isn’t actually good for performance, maybe I got the wrong idea though

There are some nuances concerning that specific statement. But the take away is: Don’t use global variables!

5 Likes

True, I just think if OPs code was written in a Julian way (rather than writing it as you’d write it in Fortran) the performance would be equivalent. I’ve never seen const + ref used in Julia (not that I’ve written anything huge in Julia I just simulate things) so I was a bit skeptical :grinning:

Yes, if the code was written more properly in Julia, the performance should be similar.

If you use const + ref to define scalars, is like defining every scalar as a one-element array, with a specific position in heap memory, carrying a pointer. Scalars do not need that, they are immutable and that allows many compiler optimizations.

3 Likes