Optimizing the Particle Life: From 400 to 4M particles

https://programmingattack.com/articles/optimizing-the-article-life-from-a-thousand-to-a-milion-particles

by caranea • 2 years ago

153 60 2 years ago

programmingattack.com

hyperbrainer3• 2 years ago

Why was C slower though?

KeplerBoy• 2 years ago

No optimization flags could be a big part of the reason.

Haven't looked closely at the code or tried it, but with -O3, -fopenmp and a well-placed pragma the performance could increase many-fold.

Heck, with NVC++ you could offload that thing to a GPU with minimal effort and have it flying at the memory bandwidth limit.

wzdd2• 2 years ago

Simply compiling it with -O3 produces something which completes in half the time of the JavaScript version (350ms for C, 750ms for JS), so perhaps that.

Edit for Twirrim: on this system (Ryzen 7, gcc 11): "-O3": 350ms; "-O3 -march=native": 208ms; "-O2": 998ms; "-O2 -march=native": 1040ms.

Edit 2: Interestingly, changing the C from float to double produces a 3.5x speedup, taking the time elapsed (with "-O3 -march=native") to 58ms, or about 12x faster than JS. This also makes what it's computing closer to the JavaScript version.

anon9463• 2 years ago

Did you add output of the results to the C code?

Twirrim4• 2 years ago

For fun and frolics:

No flags: 1843ms

-march=native: 2183 ms

-O2: 423 ms

-O2 -march=native: 250 ms

-O3: 425 ms

-O3 -march=native: 255 ms

O3 doesn't seem to be helping in my case.

jasonjmcghee2• 2 years ago

Seems crazy to me that double would produce that kind of speed up. Is float getting emulated somehow? Don't they end up the same size?

pvg2• 2 years ago

Author mentions they didn't use optimization flags but doesn't include the compilation details. You can sort of guess that (relatively) unoptimized C might perform worse than V8's JIT on short, straight computational code - you're more or less testing two native code generators doing a simple thing except one has more optimizations enabled and wins.

hyperbrainer3• 2 years ago

Oh, I assumed he still did -O2, and did not do anything else. Is that bad to assume?

PS: I do not use C beyond reading some of its code for inspiration, so kinda unaware

anon9463• 2 years ago

IMO, not using any optimization flags with C is somewhat arbitrary, since the compiler writers could have just decided that by default we'll do thing X, Y, and Z, and then you'd need to turn them off explicitly.

FWIW, without -O, with -O, and with -O4, I get 2500ms, 1500ms, and 550ms respectively. I didn't bother to look at the .S to see the code improvements. (Of course, I edited the code to output the results, otherwise, it just optimized out everything.)

hyperbrainer3• 2 years ago

(Of course, I edited the code to output the results, otherwise, it just optimized out everything.)

O(1) :)

caranea2• 2 years ago

Thanks for posting your results!

Since I was already set on writing in-browser particle life, I didn't benchmark C code with different flags.

anon9463• 2 years ago

Completely reasonable. I'd probably edit your blog post a bit to indicate that.

TylerE• 2 years ago

Should also test -Os when doing this sort of thing. Sometimes the reduced size greatly improves cache behavior, and even when not it's often outright competitive with at least -O2 anyway (usually compiles faster too!)

abainbridge3• 2 years ago

One optimization for the C code is to put "f" suffixes on the floating point constants. For example convert this line:

    t[i] += 0.02 * (float)j;

to:

    t[i] += 0.02f * (float)j;

I believe this helps because 0.02 is a double and doing double * float and then converting the result to float can produce a different answer to just doing float * float. The compiler has to do the slow version because that's what you asked for.

Adding the -ffast-math switch appears to make no difference. I'm never sure what -ffast-math does exactly.

Minimal case on Godbolt:

https://godbolt.org/z/W18YsnMY5 - without the f

https://godbolt.org/z/oc1s8WKeG - with the f

a13692099932• 2 years ago

I believe this helps because 0.02 is a double and [...] can produce a different answer

In principle, not quite. The real/unavoidable(-by-the-compiler) problem is that 0.02 is a not a diadic rational (not representable exactly as some integer over a power of two). So its representation (rounded to 52 bits) as a double is a different real number than its representation (rounded to 23 bits) as a float. (This is the same problem as rounding pi or e to a double/float, but people tend to forget that it applies to all diadic irrationals, not just regular irrationals.)

If, instead of `0.02f` you replaced `0.02` with `(double)0.02f` or `0.015625`, the optimization should in theory still apply (although missed optimization complier bugs are of course possible).

sllabres• 2 years ago

"I'm never sure what -ffast-math does exactly."

Me too when I am away from C for a while. The topic has been on HN [3]

* Enable the use of SIMD instructions

* alter the behavior regarding NaN (you can't even check for NaN afterwards with isnan(f))

* alter the associativity of expression a+(b+c) might become (a+b)+c which seems inconspicuous at first, but there are exceptions (as example see [1] under -fassociative-math)

* change subnormals to zero (even if your program isn't compiled with this option, but a library you link to your program).

A nice overview from which I summarize is in [1] which contains a link to [2] with this nice text:

"If a sufficiently advanced compiler is indistinguishable from an adversary, then giving the compiler access to -ffast-math is gifting that enemy nukes. That doesn’t mean you can’t use it! You just have to test enough to gain confidence that no bombs go off with your compiler on your system"

[1] https://simonbyrne.github.io/notes/fastmath/

[2] https://discourse.julialang.org/t/when-if-a-b-x-1-a-b-divide...

[3] https://news.ycombinator.com/item?id=29201473 (107)

CyberDildonics4• 2 years ago

This is fairly trivial stuff with underwhelming results. 4 million particles isn't much on modern hardware, even with spatial partitioning. The hardest part of spatial partitioning is probably the one dimensional partition algorithm in the first place, which can be found here:

https://en.cppreference.com/w/cpp/algorithm/partition

iainmerrick• 2 years ago

It's not rocket science, but "trivial" is very harsh. This stuff is fiddly to get right. The author (who I assume is a student) seems to have done a good job and it's a nice writeup.

CyberDildonics4• 2 years ago

I guess you could argue that, but it's strange that people want to give a student's first project (with a lot of huge mistakes like thinking node is going to outperform C because they didn't turn optimizations on) so much interest. I think people are taken in by big numbers and assume there is something cutting edge because they don't know any better.

The results are just some particles in an extremely basic pattern, there isn't a lot of payoff.

atulvi• 2 years ago

4M particles is not the question. what's the framerate?

bee_rider2• 2 years ago

You can go to the “next” story at the bottom of the post, and then the author has some links to an online implementation.

jasonjmcghee2• 2 years ago

Would vary significantly based on the GPU, I would assume.

Jupe3• 2 years ago

This is pretty neat...

Having an interest in the space for many years (on-and-off; life gets busy sometimes) - I've often lamented the lack of good pixel-level performance optimizations in graphics cards. Everything seems hell-bent on polygons.

Some years ago, DirectDraw on Windows was an excellent way to optimize the graphics portion of these types of things, but that all went away as "3D Is The Way" mentality took over.

My explorations of Processing.* were neat, but lacked the IDE and library support I like as a developer.

pixelpoet3• 2 years ago

There's no lack of "pixel-level performance optimisations" on GPUs, and just because basically the whole world seems to think that graphics programming == using some rasterisation library like OpenGL or DX or Vulkan or whatever for realtime applications, doesn't mean you're forced to use it or that that's all there is.

Just fire up OpenCL (or CUDA, if you must) and start executing some data-parallel kernels.

Jupe3• 2 years ago

That's all a little beyond me at this point. But very compelling!

As the original simulation shows, the idea is to calculate the difference and proximity of the pixels. So the drawing part is just a projection...

The ease of getting a memory buffer pointer, setting pixels within my calculation loop with simple offsets was very simple and compelling. I think I should look deeper into OpenCL/CUDA for this use-case, as I mentioned in the other response to my comment.

Thanks!

dartos• 2 years ago

Is vulkan strictly a rasterization library?

I thought you could also do numeric calculations with compute shaders.

smlacy• 2 years ago

The thing you describe as "per-pixel optimizations" is exactly what a fragment shader is? You don't need "a polygon" per say, you can just run said code over ever pixel in the frame for each frame. This is how simple screenspace effects are implemented, compositing, HDR, etc.

Jupe3• 2 years ago

Interesting. While I've played with ShaderToy a bit, I've not gotten such a thing to run locally. Maybe I should revisit this..

I guess all that I'd need would be in-memory buffer of my field, then a single function to copy the buffer (or section of the buffer if I want larger than screen fields) to the video buffer.

bee_rider2• 2 years ago

Since I’m on my phone, I can only try the 2D js versions, but the quadtree version does seem to provide some really obvious and nice frame rate improvements. It seems like good work, although I don’t know much about this field.

Sidenote:

I’m often too self conscious to put anything technical online with my real name on it, probably a bit of imposter syndrome or something. I think based on

For one of my programming classes this semester, my task was writing a program in Python. It was only an introductory course, so it didn’t have to be anything fancy.

the author is sort of early on in their career. I often suspect I’d be better off if I’d gotten over myself and thrown more early experiments online. I dunno, maybe this stuff just isn’t a big deal for some people, but I have to applaud the bravery.

_neil• 2 years ago

Give yourself to the embarrassment. Publish for your past self. Or publish for someone smarter to come along and correct you.

(Currently trying to take this advice as well)

alook• 2 years ago

Love this advice. I’ve always struggled with the same thing.

One friend started telling me to get over myself and publish, then he reworded it to:

“Get over being under yourself.”

For some reason that one stuck with me and gave me the strength to ship my first few pieces of writing.

peter_l_downs• 2 years ago

I was having trouble finding the link to the actual project, so here you go:

- Github repo https://github.com/Caranea/gpulife

- Finished particle simulator https://webgpulife.netlify.app/

danielam• 2 years ago

This reminds me of a paper[0] I co-wrote during a college internship that involved simulating the evolution of a Brownian system of particles with production and annihilation rules as the evolution of probability density distributions for each particle type within a partitioned cubic periodic box. Because we were interested in how these particles segregated spatially, visualization was a matter of finding the boundaries separating each domain of particles. I made some fun little videos of the evolution of the system over time (which tracks very conspicuously the entropy production, plotted in the paper), but unfortunately, they've since been taken down from the institute's website where they were hosted. (The paper does contain some snapshots that went into these visualizations, however.)

[0] https://pubs.aip.org/aip/jcp/article-abstract/122/17/174105/...

Twirrim4• 2 years ago

The thing that strikes me as particularly weird about this is the use of numpy there. In other languages, they're using native code. In python they're reaching out to numpy, which is a great library, but not awesome inside a hot loop unless you're keeping the operation you're carrying out within numpy itself. This means, right in that hot loop, they're doing a lot of translating of numbers between python representation and native c representation.

Cutting numpy out by making the line "t = [0.0]*l" gets it down to 17059 ms, without attempting any other optimisations. Using that plus pypy (so you get the JIT that you have with javascript/v8) gets it down to 958 ms.

To show the cost of those translations, sticking with numpy. If we switch that inner loop to:

        temp_t_i = t[i]
        temp_t_i += 0.02 * j
        temp_t_i *= 0.03 * j
        temp_t_i -= 0.04 * j
        temp_t_i /= 0.05 * (j+1)
        temp_t_i = temp_t_i

It speeds us up from 47552 ms to 28251 ms. Almost half the execution time. That's still doing two hops back and forth, though. If you cut it down to a single line it's even faster at 18458 ms, cutting execution time down to about a third of the original example. Pypy isn't able to help here at all, this is sort of a pathological case for it.

edit: I'll add, I'm not that good with numpy, rarely use it myself. Not sure if it's possible to do that inner loop all within numpy somehow. I imagine that'd be a lot faster still.

nerpderp82• 2 years ago

This is the problem that Numba was designed to solve.

https://numba.readthedocs.io/en/stable/user/5minguide.html

cl3misch2• 2 years ago

No. The problem in this post can be vectorized (== expressed as array OPs) with idiomatic numpy, making it very fast (see sibling comment).

Numba is designed to speed up code where looping is unavoidable, i.e. the code can't be (easily) expressed as array OPs.

bscottmay• 2 years ago

Yeah, an all-numpy version runs in less than 1ms on my M1 air

  import numpy as np
  l = 10_000
  t = np.empty(l, dtype=np.float32)
  j = np.arange(l)
  t = 0.02 \* j
  t *= (0.03 * j)
  t -= (0.04 \* j)
  t /= 0.05 \* (j + 1)

Twirrim4• 2 years ago

So rather than spend 10-20 minutes reading about numpy, the author wrote 3 other implementations...?

The fact that they ran the C code without the optimisation flags and compared it that way makes me think Javascript was what they actually wanted to write this one in anyway.

chris_va• 2 years ago

Sounds like they are slowly learning electrostatic PIC codes from first principles.

Optimizing the Particle Life: From 400 to 4M particles

Related Stories